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ATMOSPHERIC  TURBULENCE  MODELLING  IMPACTING  STRATOSPHERIC  ENVIRONMENT  PROGRAM 


OBJECTIVES : 

The  primary  objectives  of  this  study  can  be  described  as  follows: 

1)  Develop  methods  of  describing  the  constituents  and  thermal 
distributions  of  the  normal  stratosphere  and  calculate  how  these 
distributions  are  affected  by  Air  Force  flight  operations  and 
atmospheric  transport  coefficients. 

2)  Characterize  stratospheric  disturbances  that  can  be  produced 
from  mountain  waves  and  underlying  thunderstorm  activity. 

3)  Develop  specialized  numerical  techniques  to  Incorporate  radiation 
transfer  and  couple  photochemical  processes  with  the  heat  balance 
of  the  stratosphere. 

This  task  Is  the  key  to  the  stratospheric  environment  project  wherein 
the  results  of  the  various  measurements  are  used  to  establish  a  realistic 
model  of  the  stratosphere.  The  model.  In  turn.  Is  used  to  predict  If 
environmental  changes  will  occur  as  a  result  of  AF  flight  operation.  This 
task  is  essential  to  the  preparation  of  environmental  impact  statements 
for  the  Air  Force. 

The  major  effort  of  this  work  on  turbulence  modelling  Is  to  Improve 
stratospheric  models  for  the  prediction  of  Air  Force  systems  operations  in 
the  stratosphere.  In  particular,  the  probability  of  occurrence  of 
turbulence  for  a  given  location  and  season  as  a  function  of  altitude  is  a 
required  Input  to  the  tropospheric/stratospheric  model  used  to  predict 
effects  of  Air  Force  Systems  emission  pollutants  in  the  stratosphere.  Thus 
we  are  investigating  the  seasonal  variability  of  the  characteristic  layers 
of  turbulence  In  the  troposphere  and  stratosphere  and  determining  the 
seasonal  variations  in  the  natural  levels  of  NO^  trace  constituents  which 
are  also  engine  combustion  by-products.  This  will  result  In  refinements 
in  the  troposphere/stratosphere  1-D  and  3-D  chemical/radiation  transport 
models  which  are  used  to  predict  environmental  effects,  such  as  ozone 
depletion,  caused  by  Air  Force  Systems  emissions  in  the  troposphere  and 
stratosphere. 
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Meterologlcal  data  (temperature,  winds,  pressure  and  relative 
humidity)  in  the  height  range  0-60  km  on  magnetic  tape  is  currently  being 
restructured  to  provide  efficient,  easy  access  for  processing  with 
minimum  computer  time.  A  statistical  approach  is  being  used  to  determine 
the  seasonal  probability  of  occurrence  of  turbulence  as  a  function  of 
latitude,  longitude  and  altitude  through  the  use  of  the  Richardson  criteria 
for  the  presence  or  absence  of  turbulence  in  the  atmosphere.  NO^  trace 
constitutents  measured  by  LKD  will  be  compared  with  data  from  other 
experimental  techniques  to  determine  the  extreme  of  seasonal  variability 
in  these  constituents.  The  turbulence  and  NO  seasonal  results  will  be 

X 

used  to  parameterize  the  Important  role  of  transport  processes  and  provide 
fine  adjustments  in  the  tropospheric/stratospheric  chemistry  of  the  1-D 
and  3-D  models . 


I 
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INTRODUCTION : 


The  stratospheric  seasonal  variations  task  under  the  Stratospheric 
Environment  Project  is  currently  studying  the  general  aspects  of  turbulence 
in  the  troposphere  and  lower  stratosphere.  The  data  being  used  in  this 
study  has  been  provided  by  the  Environmental  Data  Service  National  Climatic 
Center  in  Asheville,  No.  Carolina.  The  "rawinsonde  data"  consists  of 
winds,  temperature  and  pressure  as  a  function  of  altitude  obtained  from  all 
available  stations  (currently  144  stations)  for  the  period  1948-1976.  At 
the  present  time  only  the  1970-1976  data  have  been  reformated  for  process¬ 
ing  on  the  CDC  6600  computer  system.  This  data  is  being  used  in  calculating 
the  Richardson  number,  ^  stability  criterion  for  determining  the  presence 
or  absence  of  atmospheric  turbulence,  and  is  given  by  the  following 
relationship: 

'1-*^  [-|f^ 

g  is  the  acceleration  of  gravity 

6  T/ c  is  the  vertical  temperature  gradient  . 
o  z 

r  is  the  dry  adiabatic  lapse  rate  (  10°  ) 

[""It]  '^Vy/gp2j^s  the  square 

of  the  vertical  shear  of  the  horizontal  winds. 


where:  " 


and 
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The  rawinsonde  system  is  not  ideally  suited  to  take  measurements  to 
provide  thermodynamic  data  at  the  required  sampling  intervals  for  this 
application.  The  data  being  used  in  this  study  was  originally  intended  to 
provide  northern  hemisphere  pressure  maps  for  a  number  of  millibar  pressure 
levels.  In  order  to  determine  the  Richardson  number  at  equal  height 
Intervals,  the  temperature  and  wind  component  data  points  are  fit  using  a 
Hermite  Interpolation  algorithm.  The  temperature  function  along  with  the 
derivatives  of  the  winds  and  temperature  with  respect  to  altitude  are 
Interpolated  at  1/10  km  Intervals  and  used  to  calculate  the  Richardson 
number.  The  probability  of  occurrence  of  turbulence  is  obtained  by 
defining  a  critical  Richardson  number  (R^  *  1  to  1/A  according  to  current 
literature)  and  assuming  turbulence  present  when  the  Richardson  number 
is  equal  to  or  less  than  this  value.  The  probability  is  then  taken  as  the 
ratio  of  the  total  number  of  occurrences  to  the  total  number  of  measure¬ 
ments  at  a  1.0  km  level  over  a  season.  These  seasonal  results  are  then 
averaged  over  one  km  and  plotted  as  in  the  figure  provided.  Two  techniques 
are  used  to  establish  the  seasonal  tropopause  height.  The  seasonally 
averaged  temperature  gradient  value  is  tested  for  a  sustained  change  in 
sign  which  indicates  the  change  in  lapse  rate  from  the  troposphere  to  the 
stratosphere.  Frequently,  at  high  latitudes,  the  temperature  remains 
fairly  constant  above  the  troposphere  and  through  the  stratosphere.  Since 
a  negative  lapse  rate  may  not  develop,  a  second  test  is  used  which 
established  the  tropopause  height  as  the  maximum  rate  of  change  in  the 
temperature  gradient. 

Initial  results  from  stations  at  low,  middle  and  high  latitudes 
reveal  interesting  seasonal  and  latitudinal  variation  in  the  occurrence  of 
turbulence  as  well  as  particular  altitude  levels  which  consistently  reveal 
relatively  high  occurrences  of  turbulence.  It  is  also  interesting  to  note 
that  turbulence  rapidly  decreases  at  heights  which  generally  agree  with 
those  levels  reported  in  models  as  the  tropopause  level  and  follows  this 
pattern  as  the  tropopause  level  decreases  from  its  highest  level  in 
equatorial  regions  to  its  lowest  levels  at  high  latitudes.  This  particular 
result  is  substantiated  in  theory  since  the  stratosphere  region  above  the 


tropopause  level  Is  a  dry  stable  region  having  a  negative  temperature 
lapse  rate,  and  thus  the  degree  of  turbulent  mixing  and  transport  Is 
greatly  reduced  from  that  occurring  in  the  troposphere. 

The  resolution  of  the  temperature  and  wind  data  (100m  to  a  few 
kilometers)  varies  with  the  rise  rate  of  the  rawinsonde  balloon  and  the 
criteria  of  selecting  data  points  between  mandatory  levels  (those  standard 
millibar  height  levels  recorded  at  all  stations).  The  fit  applied  to  the 
data  to  provide  Richardson  numbers  at  equal  altitude  Intervals  brings  Into 
question  the  validity  of  the  Richardson  number  determination  by  possibly 
Introducing  effects  not  present  in  the  thermodynamic  structure  of  the 
atmosphere.  The  possibility  of  obtaining  a  limited  amount  of  data  at 
suitable  sampling  intervals  to  obtain  some  knowledge  of  the  actual 
structure  in  the  data  has  been  Investigated.  Personnel  at  the  Chatham,  MA. 
rawinsonde  station  have  cooperated  In  the  Initial  Investigation  by 
demonstrating  the  present  techniques  in  data  acquisition  during  a  visit 
to  chat  site  in  September,  1978.  They  also  provided  at  that  time  an  analog 
strip  chart  recording  of  temperature  along  with  the  standard  printout  of 
thermodynamic  data.  Data  points  extracted  from  the  temperature  analog  at 
the  smallest  possible  interval  indicate  spectrum  wave  lengths  of  approxi¬ 
mately  a  kilometer  or  longer. 
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TURBULENCE  STimv  SOFTWARE 


A  description  of  the  software  used  In  this  study  is  presented  next 

This  is  followed  with  a  summary  of  statistical  procedure  which  we  have  de¬ 
veloped  . 
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turbulence  study  flowchart 


FOUKSEASON  PLOT 


The  primary  function  of  this  plot  is  to  use  the  Season  E  output  and 
plot  the  Percent  Occurrence  of  Turbulence  (both  RI  £  .25  and  RI  £l.O) 
versus  altitude.  Each  plot  contains  one  year’s  data  with  all  four  seasons, 
and  altitude  2  to  25  kilometers.  Each  plot  presents  the  station  name, 
number,  longitude,  latitude  and  year.  Also  incorporated  into  this  program 
is  the  determination  of  the  seasonal  tropopause  levels.  This  is  accomplish 
ed  and  displayed  using  two  approaches.  Method  one  tests  for  a  sustained 
sign  change  in  the  temperature  gradient.  Method  two  looks  for  the  maximum 
change  in  slope  of  the  temperature  gradient.  Output  is  on  both  pen  plots 
and  microfiche. 


LATPLOT 


This  program  uses  data  extracted  from  the  Season  E  output.  Created 
is  a  plot  of  achosen  parameter  versus  latitude.  Numerous  stations  arranged 
by  latitude  can  be  shown  on  one  plot.  The  desired  parameter  with  its 
standard  deviation  is  plotted  for  a  chosen  altitude.  Any  or  all  seasons 
can  be  included  on  one  plot.  This  plot  is  very  good  in  showing  relation¬ 
ships  between  station,  and  through  latitudinal  changes. 

DATABASE 


The  data  base  program  examines  the  Percent  Occurrence  of  Richardson 
Number  greater  than  -4.0  and  less  than  or  equal  to  1.0.  This  program 
could  also  process  any  of  the  other  parameters.  Our  analysis  is  done 
seasonally  for  years  1971  through  1975  and  kilometer  altitude  levels  2-25. 
Calculations  are  performed  by  two  methods,  a  five  year  average  for  each 
kilometer  level  and  an  average  over  the  total  2-25  kilometer  range  per  year 
Values  calculated  are  the  mean,  standard  deviation,  maximum,  minimum  and 
range.  We  also  perform  a  two-way  analysis  of  varirance.  Results  include 
the  overall  mean,  variance  component  due  to  altitude  bins,  variance 
components  due  to  yearly  variations,  total  variance  component,  residual 
error  and  standard  error.  These  are  done  for  the  entire  2-25  km  range  and 
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also  for  altitude  groups  2-7,  8-13,  14-19  and  20-25. 

Output  consists  of  for  each  season,  the  24  levels  of  yearly  means  and 
the  24  levels  of  the  5  year  means  on  permanent  file  storage. 

BASEPLOTX 


This  plot  uses  the  yearly  data  base  output  plus  the  5  year  average 
output.  It  created  plots  of  Percent  Occurrence  of  Turbulence  (RI  LE.l) 
versus  altitude.  A  set  of  plots  is  created  for  each  season  with  each  of  the 
5  years  and  the  5  year  average  plotted  per  set.  The  sets  are  small  and 
concise,  ideal  for  comparing  yearly  and  seasonal  variations.  Both  pen  and 
microfiche  plots  are  available. 

LOWALT 


The  Lowalt  program  selects,  through  card  input,  the  station  and  year 
desired.  The  program  searches  the  OS  tape  and  extracts  the  raw  data.  The 
extracted  records  consist  of  a  leader  containing  deck  number,  station 
number,  year,  month,  day,  hour  and  number  of  levels  recorded.  Each  level 
contains  the  following  parameters  -  pressure,  altitude,  temperature, 
relative  humidity,  wind  speed  and  direction  and  level  number.  The  data 
has  an  altitude  range  limit  of  35  kilometers.  Output  is  a  7  track  CC 
tape  containing  the  information.  Approximately  15  years  data  recorded 
per  tape. 


ORIGINAL  OS  TAPES 


The  original  OS  tapes  are  a  set  of  40  9-track  binary  tapes.  They 
contain  Rawlnsonde  data  from  144  stations  throughout  the  Western  Hemisphere. 
Testing  occurred  from  1970  to  1976. 
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ALTFIND 


Altflnd  is  a  program  created  to  examine  the  altitude  structure  of 
the  Season  E  input.  We  take  the  validated  Input  points  from  each  launch 
and  segment  It  Into  altitude  levels  of  2-7.5,  7.5-13.0,  13.0-20.0  and 
20.0-25.0  km.  We  find  the  number  of  points  and  the  average  distance 
between  them  for  each  range.  We  then  take  a  seasonal  average  of  these 
values . 


SEASON  E 


Season  E  reads  the  Lowalt  generated  7  track  CC  tape.  The  data  for 
each  launch  (2  per  day)  Is  verified  and  processed  at  each  whole  kilometer 
level  by  a  spline  interpolation.  The  desired  turbulence  parameters  are 
calculated  and  averaged  by  season.  The  variables  of  major  Interest  are 
the  Richardson  number  and  the  percent  occurrences  of  Richardson  number 
less  than  .25  and  less  than  1.0.  Output  is  a  permanent  file  containing 
the  38  variables  at  each  kilometer  level,  1  to  35,  for  all  four  seasons 
and  years  1970-1976.  Another  output  option  presents  the  daily  raw  input 
values  and  the  dally  calculated  turbulence  parameters. 

BASEPL0T2PL 

This  plot  uses  the  5  year  average  data  base  output.  Two  plots  for 
each  season  are  created.  Plot  one  shows  the  5  year  average  per  altitude 
data.  The  mean  Is  plotted  and  the  standard  deviation,  maximum  and  minimum 
are  shown  on  the  side.  The  second  plot  shows  the  yearly  means  taken  over 
the  total  altitude  range.  The  plot  is  good  to  examine  yearly  variations. 
The  yearly  standard  deviation,  maximum  and  minimum,  are  also  shown. 

The  main  parameter  plotted  has  been  the  Percent  Occurrence  of 
Turbulence  although  many  of  the  other  variables  have  been  examined  through 
these  plots.  Output  pen  plots  or  microfiche. 
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RJBGRESSION  PACKAGE 


The  Regression  Package  is  an  extremely  important  element  in  the 

analysis.  We  use  the  data  base  output  as  input,  and  through  this  package 

perform  a  wide  variety  of  statistical  analysis.  We  examine  the  percent 

2  3 

occurrence  with  parameters  altitude,  altitude  ,  altitude  ,  latitude, 

2  2 

latitude  ,  longitude,  longitude  ,  the  four  seasons  and  the  five  years. 
With  this  package  we  can  examine  the  effects  of  each  of  these  parameters 
on  our  model.  Output  includes  means,  deviation,  correlation  coefficients, 
corrected  sums  of  squares  and  cross  products,  the  sum  of  the  residuals 
squared,  multiple  correlation  coefficient,  standard  error  of  estimates, 
degrees  of  freedom,  F-values  and  more.  Scatter  plots  are  also  created. 

Typical  percent  occurrence  plots  follow. 
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TURBULENCE  STUDY  -  STATISTICAL  ANALYSIS 


OBJECTIVES: 

The  statistical  analysis  plan  to  be  described  below  has  two  major 

obje-'.;  Live'i  : 

1)  To  relate  occurrence  of  turbulence  to  latitude,  longitude, 
altitude  and  seasonal  variation, 

2)  To  identify  areas  of  high  turbulence. 

I  -  SITE  SELECTION 

At  present  data  are  available  on  1A4  sites  (locations)  for  the  years 
1971-1976,  From  these  sites  we  have  selected  a  subset  of  14  sites  for 
initial  analysis.  This  subset  has  been  selected  to  include: 

1)  A  wide  range  of  latitude  and  longitude.  The  14  sets  are  repre¬ 
sentative  of  sites  ranging  from  70  to  120  longitude  (W)  and 
25  to  50  latitude  (N) .  This  will  give  the  results  broad 
applicability. 

2)  "Possibly”  interesting  variables  such  as  land  configurations  (i.e. 
coastal,  mountains,  flat  lands).  This  will  allow  us  to  investi¬ 
gate  variables  that  may  introduce  variability  into  our  results. 

3)  Stations  which  are  geographically  close.  This  will  allow  us  to 
investigate  the  correlation  between  nearby  sites. 

4)  Stations  that  also  have  data  from  1948  to  1970.  This  will  permit 
a  more  extensive  time  series  analysis  (if  deemed  useful). 

The  results  obtained  from  the  14  selected  sites  should  be  representa¬ 
tive  of  locations  in  the  continental  United  States.  When  the  analysis  is 
completed  on  these  sites  data  from  other  sites  will  be  analyzed  to  act  as 
a  validation  of  the  results  obtained  from  the  14  original  sites. 

II  -  DATA  FOR  ANALYSIS 

Dependent  Variable  -  The  major  dependent  variable  is  the  proportion 
of  days  for  a  given  time  period  (e.g.,  a  season)  with  turbulence.  Turbu¬ 
lence  is  determined  according  to  whether  the  Richardson  number  is  less  than 
.25.  At  each  site  we  have  for  a  given  balloon  launch  at  most  one  Richard¬ 
son  value  for  each  1.0  k  interval  from  ground  level  up  to  a  height  of  35km, 
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A  typical  situation  is  to  have  approximately  two  balloon  launches  per  day 
providing  two  readings  for  up  to  25  km.  At  this  point  the  balloon  usually 
bursts  and  measurements  are  no  longer  available. 


Independent  Variable  -  The  independent  variables  Include:  site 
latitude,  site  longitude,  altitude,  average  temperature,  north-south  wind 
speed,  east-west  wind  speed,  temperature  variation  (e.g,,  range  of  tempera¬ 
ture  over  a  season),  density  humidity,  year,  and  season  (i.e.,  appropriate 
time  period) . 

Ill  -  ANALYSIS  OF  INDIVIDUAL  SITES  (FIRST  ANALYSIS) 

The  first  analysis  is  a  site-by-site  analysis.  That  is,  each  site 
is  analyzed  Individually. 


A  -  Summary  Statistics 

For  a  given  time  period  (season)  summary  statistics  on  all 
variables  will  be  computed.  These  include  means ,  standard 
deviations,  maximums/minimums  and  ranges  for  each  altitude 
bin  (each  1  kilometer  is  referred  to  as  a  bin) . 

In  particular,  we  determined  the  number  of  altitude  bins 
that  can  be  analyzed  with  good  statistical  accuracy  and  pre¬ 
cision.  High  altitude  bins  (say  above  25  km)  often  did  not 
have  sufficient  data  for  a  good  analysis.  Based  on  the 
results  we  obtained  from  the  14  selected  sites,  we  determined 
how  many  bins  should  be  carried  further  in  our  analysis. 

B  -  Altitude  Analysis 

For  each  altitude  bln  and  time  period  (season)  a  two  way 
analysis  of  variance  will  be  performed  using  Time  Period  by 
year  as  the  variables  for  classification.  The  objectives  will 
be:  (1)  quantify  yearly  variation  and  (2)  quantify  time 

period  variation. 

C  -  Altitude  by  Time  Period  by  Years  Analysis 

This  analysis  quantified  variations  due  to  altitude,  variations 
due  to  years,  variations  due  to  time  period  and  Interactions 
among  these  three  sources. 
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D  -  Regression  Analysis 

The  regression  model  related  proportion  of  turbulence  to  year, 
time  period,  altitude  and  weather  variables  such  as  temperature, 
density,  humidity,  wind  direction,  and  Interaction  of  these 
variables.  Specifically  the  Independent  variables  are:  Year, 
Time  Period  (season  or  month).  Altitude,  Average  Temperature, 
Temperature  Variation,  Wind  Direction  (North-South,  East-West), 
Density,  Humidity  and  Interaction  of  these. 


IV  -  COMBINATION  OF  SITES  (SECOND  ANALYSIS) 

The  second  analysis  consists  of  combining  sites  and  analyzing  them 
jointly.  If  deemed  necessary  the  following  analysis  can  be  performed  on 
subsets  of  sites. 

A  -  Regression  Analysis  (altitude  x  latitude  x  longitude) 

A  regression  relating  proportion  of  turbulence  to  geographic 
location  and  altitude  will  be  performed.  This  regression  will 
also  include  year  and  time  period  (season) .  The  model  is 

(Proportion  turbulence)  =  b^  +  b^  (Altitude)  +  b2  (latitude) 

+  b^  (longitude)  +  b^  (year)  +  b^  (time  period). 

If  necessary  the  model  will  be  expended  to  include  interactions 
of  the  above  variables. 

B  -  Expanded  Regression  Analysis 

The  next  analysis  will  further  expand  on  the  above  regressions. 
Variables  which  quantify  land  configurations  will  be  added. 

Also  weather  variables  as  discussed  in  III  D  will  also  be  added. 

V  -  LOGISTIC  REGRESSION 

The  objective  of  analysis  is  to  relate  occurrence  of  turbulence  to 
sets  of  variables.  In  particular  we  want  to  relate  variables  to  the 
probability  of  turbulence.  A  useful  model  for  this  Is 
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'  ^  e  b^Xi+b 
1  +  e  1 


o 


where  P  represents  the  probability  of  turbulence  and  E  b^^Xi+b^ 

represents  a  linear  function  relating  variables  Xj^,...,X^  to  the  probability 

of  P. 
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By  use  of  the  tranafonuation 


■■  w 

we  have  In  (  - ^  ^  )  =  b  +  E  h  X 

\  W  oil 

and  regression  analysis  is  now  possible  using  the  transformed  dependent 
variable.  We  will  perform  these  logistic  regressions  in  addition  to  the 
above  regressions. 


A  -  Individual  Sites 

A  logistic  regression  will  be  performed  using  variables  listed 
in  III  D. 

B  -  Combination  of  Sites 

Logistic  regressions  will  be  performed  using  variables  described 
in  IV  A  and  also  IV  B. 

Note:  Analysis  IV  A,  IV  B,  V  B  will  allow  us  to  relate  the  occurrence  of 
turbulence  to  latitude,  longitude,  altitude  and  seasonal  variation. 
Additionally  these  analyses  will  quantify  Importance  of  weather 
variables  to  turbulence.  The  functions  developed  in  these  analyses, 
especially  the  results  of  V  B,  will  permit  us  to  identify  the  profiles 
of  variables  that  relate  to  higher  turbulence  -  viz,  we  will  be  able 
to  identify  values  of  X^,...,X  given  in  V  that  produce  high 
probabilities  P  given  by  the  logistic  regression. 

Note:  While  the  above  analyses  address  mainly  the  problem  of  explaining 
past  data,  we  will  also  investigate  their  ability  to  predict 
turbulence . 

Note:  As  stated  earlier  the  above  will  first  be  performed  on  the  14  selected 
sites.  When  this  is  completed  other  sites  in  the  continental  U.S. 
will  be  analyzed  for  validation.  Then  high  latitude  sites  will  be 
selected  for  an  analysis  similar  to  the  one  described  above. 
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Rocketsonde  Measurements  -  Radiosonde  observations,  or  "RAOBS"  as  they  are 
termed,  are  made  to  determine  the  pressure,  temperature,  and  humidity  from 
the  surface  to  the  point  where  the  balloon  bursts.  The  radiosonde  consists 
of  meteorological  measuring  elements  coupled  to  a  radio  transmitter  and 
assembled  into  a  small  lightweight  box.  The  device  is  carried  aloft  by 
a  balloon  filled  with  hydrogren  or  helium  gas.  Included  in  the  train  is 
a  small  parachute  to  slow  the  descent  of  the  instrument  after  the  balloon 
bursts,  thereby  minimizing  the  danger  of  injury  to  life  and  property. 

As  the  balloon  rises,  measurements  of  pressure,  temperature,  and  relative 
humidity  are  transmitted  to  a  gound  station  where  the  data  are  recorded 
automatically.  An  observer  then  transcribes  the  information  into  a  more 
commonly  used  form  acnd  plots  it  on  various  charts.  Measurements  of 
pressure  are  made  in  millibars,  temperature  in  degrees  Celsius,  and  moisture 
in  percent  of  relative  humidity.  The  Celsius  temperature  scale  is  the 
same  as  the  centigrade  scale  and  replaces  that  scale  in  all  meteorlogical 
observations. 

As  the  balloon  rises,  it  is  followed  either  visually  by  an  observer 
at  a  theodolite,  or  electronically  by  radlo-dlrection-flnding  equipment 
that  tracks  the  transmitted  signal.  The  balloon's  drift  away  from  the 
release  point  is  plotted,  and  from  this  the  direction  and  speed  of  the 
air  movements  are  determined.  The  wlnds-aloft  observation  is  termed 
a  "RABAL"  when  the  tracking  is  done  visually,  and  a  "RAWIN"  when  the 
tracking  is  done  electronically.  A  combined  rawln  and  raob  is  termed  a 
"RAWINSONDE". 


Rawlnsonde  System  -  The  basic  rawinsonde  system  consists  of  a  balloon- 
borne  radiosonde,  a  receiver  and  tracking  unit,  and  a  recorder.  The  most 
commonly  used  radiosonde  transmits  meteorlogical  information  consisting  of 
pressure  obtained  from  an  aneroid  cell,  temperature,  and  relative  humidity. 
The  telemetered  meteorological  information  from  the  tadiosonde  is  detected, 
amplified,  and  shaped  by  a  receiver  and  the  processed  information  is  print¬ 
ed  in  graphic  form  on  a  strip  chart  recorder.  An  antenna  array,  which  is 
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an  integral  part  of  the  receiving  equipment,  is  used  to  track  the 
radiosonde.  Tracking  is  accomplished  either  manually  or  automatically. 

The  elevation  and  azimuth  angle  information  so  obtained,  in  conjunction 
with  height  information  computed  from  the  meteorological  data,  is  used  to 
compute  winds  aloft. 

Radar  Wind  -  A  target,  suitably  mathced  to  the  frequency  of  the  radar  in 
use,  is  attached  to  the  train  of  the  radiosonde  baloon.  The  radar  will 
provide  the  slant  range  to  the  target  by  measurement  of  the  time  interval 
between  the  transmitted  and  reflected  radar  pulse.  The  azimuth  Indicator 
of  the  radar  set  will  indicate  the  azimuth  angle  of  the  target.  Slant 
range  is  used  in  conjunction  with  height  data,  computed  from  the  radio¬ 
sondes'  meteorological  data,  to  compute  distance  from  observation  point. 

Radiosonde  Instruments  -  The  radiosonde  is  a  baloon-borne ,  battery-powered 
instrument  used  together  with  the  ground  receiving  equipment  to  delineate 
the  vertical  profile  of  the  atmosphere.  Pressure  is  measured  by  means  of 
a  baroswitch  which  employs  an  expanding  aneroid  pressure  cell  to  move  a 
contact  arm  across  a  commutator  bar  as  the  pressure  decreases.  Temperature 
is  measured  by  a  thermistor,  the  electrical  resistance  of  the  thermistor 
being  a  function  of  temperature.  Relative  humidity  is  measured  by  a 
hygrlstor,  the  electrical  resistance  of  the  hygrlstor  being  a  function  of 
relative  humidity  and,  to  some  extent,  temperature.  As  the  radiosonde 
ascends,  the  thermistor  and  hygrlstor  are  switched  sequentially  into  the 
modulator  circuit  by  the  baroswitch.  The  frequency  of  the  received  signal, 
therefore,  is  alternately  a  function  of  temperature  and  humidity,  and  may 
be  any  value  from  10  to  200  Hz.  In  addition,  fixed  resistors  are  connected 
periodically  into  the  modulator  circuit  by  the  baroswitch  so  that  the 
frequency  transmitted,  when  these  are  connected  will  show  any  change  that 
may  be  occurring  In  the  modulator  characterist ic  due  to  variations  of 
circuit  parameters.  Radiosondes  are  provided  with  the  thermistor  mounted 
externally  on  an  outrigger  which  exposes  it  at  some  distance  from  the 
radiosonde  case. 
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The  hygrietor  Is  enclosed  within  an  alrduct  to  shield  the  element  from 
direct  exposure  to  precipitation.  The  radiosonde  contains  an  aneroid 
baroswitch  for  measuring  pressure,  a  white-coated  thermistor  for  tempera¬ 
ture,  and  a  carbon  hygristor  for  relative  humidity.  Radiosondes  of  this 
type  will  be  considered  the  basic  radiosonde  and  will  usually  be  identified 
by  a  serial  number  containing  no  suffix. 
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Defense  Meteorological  Satellite  la  equipped  with  a  swept-f requency 
HF  noise  receiver.  The  receiver  provides  measurements  of  radio  noise  of 
terrestlal  origin  every  100  KHz  In  the  frequency  range  of  1.2  to  13.9  Mhz. 
Satellite  is  basically  In  a  sun-synchronous  morning  or  evening  orbit.  By 
analyzing  data  from  successive  morning  or  evening  orbit,  global  maps  of 
the  noise  Intensity  can  be  obtained  at  the  satellite  altitude.  Work  Is 
being  done  In  order  to  obtain  world  maps  of  noise  Intensity  so  as  to  assess 
the  morphological  behavior  of  the  noise  at  satellite  altitude  as  a  function 
of  frequency.  Our  studies  have  been  concentrated  on  the  frequencies  of  1.5 
to  13.5  Mhz  in  steps  of  0.5  Mhz  at  both  the  dawn  and  dusk  orbits.  Monthly 
average  noise  maps  are  obtained  In  order  to  provide  more  representative 
maps  of  the  average  noise  behavior. 

The  Air  Force  Global  Weather  Central  processes  several  types  of 
special  sensor  data  including  the  data  from  the  swept-f requency  HF  noise 
receiver.  The  DMSP  tapes  have  been  processed  by  GWC  on  the  UNIVAC  1110 
system,  where  the  word  length  is  36  bits. 

We  next  present  a  description  of  the  overall  software  used  in  this 
study.  Following  this  is  a  report  deriving  the  optimal  procedure  for 
estimating  the  FoF2  of  the  Ionosphere. 
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IONOSPHERIC  RADIO  NOISE  MEASUREMENT 
SOFTWARE 
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SOFTWARE  DESCRIPTION 


PROGRAM;  RUSH 

RUSH  unpacks  data  from  one  DJ^P  satellite  tape.  UNiyAC  record  Is 
read  and  put  into  CDC  words.  All  frequencies  between  1,5  and  13.5  Mhz 
in  steps  of  0.5  Mhz  are  translated  and  loaded  Into  frequency,  latitude, 
and  longitude  matrix. 

PROGRAM;  BRUSH 

Intermediate  step  before  statistical  analysis.  Data  files  created 
by  RUSH  are  combined  until  monthly  matrix  Is  complete.  Up  to  five  data 
seto  will  be  combined  with  an  existing  data  file. 

PROGRAM;  PIS 

Program  displays  average  count  matrix  for  dawn  and  dusk  and 
statistics  about  these  matrices.  Also  printer  contour  Is  done  to  provide 
high  and  low  noise  peak  for  each  frequency,  at  each  latitude  and  longitude. 

Data  Base  generated  now  used  for  further  analysis  such  as  estimates 
foF2  as  described  in  this  report. 
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I  INTRODUCTION 


This  report  describes  an  estimation  procedure  for  determining  the 
critical  frequency  of  the  ionospheric  F2  layer.  The  estimate  is  obtained 
from  measurements  of  radio  noise  obtained  by  a  satellite  orbiting  above 
the  F2  layer.  A  spectrum  analyzer  obtains  noise  power  as  a  function  of 
frequency  in  a  band  wide  enough  to  include  all  possible  values  of  the 
critical  frequency. 

The  power  measured  at  a  particular  frequency  is  a  function  of  the 
radio  energy  being  transmitted,  and  the  transmittance  of  the  F2  layer.  The 
ionospheric  parameters  depend  on  the  position  of  the  satellite  (particularly 
latitude),  local  time,  and  solar  activity.  The  noise  measured  by  the 
satellite  receiver  is  considered  to  consist  of  two  parts.  The  first  is 
referred  to  as  terrestrial  noise,  which  passes  through  the  F2  layer. 

This  noise  originates  from  such  sources  as  radio  transmitters  located  near 
the  subsatellite  point.  The  other  component  is  background  noise,  which  is 
not  affected  by  the  F2  layer,  and  can  originate  from  within  the  satellite 
instrumentation  itself  or  other  sources  located  above  the  F2  layer. 

The  available  data  comes  from  a  DMSP  satellite  with  an  on-board 
swept-frequency  receiver.  Counts  arc  obtained  for  128  frequency  channels 
from  1.2  MHz  to  13.9  MHz  in  increments  of  100  KHz.  One  important  aspect 
of  the  data  is  that  there  is  coupling  between  adjacent  channels,  that  is, 
a  fraction  of  the  power  at  one  particular  channel  is  also  picked  up  by 
adjacent  channels.  This  is  most  likely  caused  by  the  bandwidth  of  each 
channel  filter  being  wider  than  the  Interchannel  separation.  While  this 
produces  some  smoothing  of  the  data,  it  also  modifies  the  probability 
distribution  of  the  data  from  what  one  would  normally  expect.  Thus,  it 
must  be  dealt  with  before  the  model  of  the  data  is  constructed. 
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For  estimation  of  the  critical  frequency  foF2,  the  maximum  likelihood 
(ML)  technique  Is  the  method  of  choice.  It  Is  a  well  known  and  studied 
estimation  method.  For  most  applications  It  Is  well  behaved  and  the 
optimum  method.  It  Is  also  systematic  and  Intuitively  meaningful  since  It 
can  be  derived  directly  from  the  model  of  the  data. 
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II  THEORY  OF  ESTIMATION  METHOD 


A.  The  Model 

The  effect  of  the  foF2  layer  on  the  radio  noise  spectrum  received  by 
the.  satellite  can  be  modelled  as  a  high  pass  filter.  There  are  many 
candidate  filters  which  can  be  used.  For  our  purposes,  it  is  desirable 
to  use  a  filter  with  the  least  number  of  parameters,  since  each  parameter 
must  either  be  estimated  or  assigned  a  value.  A  suitable  filter  is  the 
Butterworth  filter  of  order  n. 


»n 


1  +  (- 


fo 


2n 


(1) 


The  parameter  fo  thus  becomes  the  ionospheric  critical  frequency.  The 
filter  order  n  must  also  be  chosen.  Our  procedure  will  be  to  choose  several 
values  for  n  and  observe  Its  effect  on  the  estimation  procedure. 

The  satellite  receiver  measures  the  power  at  discrete  frequencies. 

A  record  of  count  data  can  be  expressed  as 


S(f,) 


H  (f,  ) 
n  k 


for  k  «  1, ,  . .  ,N 
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where  S|^(fj^)  is  the  noise  power  due  to  background  and  instrumentation 
noise,  and  S^(fj^)  is  the  noise  power  due  to  radio  signals  of  terrestrial 
origin  which  are  affected  by  the  F2  layer.  However,  and  are  random 
processes  whose  moments  are  unknown  functions  of  frequency.  In  the  absence 
of  such  information  about  these  signals,  we  need  to  make  certdln assumptions; 
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namely,  that  these  moments  are  constant  with  respect  to  frequency.  In 
actuality,  we  need  only  consider  one  of  the  moments  (such  as  the  mean  or 
variance)  since  these  signals  are  chi-square  distributed  with  all  moments 
dependent  on  one  of  the  moments.  Thus,  we  obtain  one  additional  parameter 
for  Sjj  and  one  for  which  must  be  estimated,  along  with  the  Important 
quantity  fo. 

The  problem  can  thus  be  expressed  as  a  classical  parameter  estimation 
problem.  The  procedure  Is  as  follows: 

1.  Preprocess  the  data  Into  a  form  acceptable  by  the  estimation 
method,  described  in  part  B. 

2.  Estimate  the  critical  frequency  fo  (and  the  and 
parameters  as  a  by-product)  by  the  maximum  likelihood  method 
(ML) ,  described  in  part  C. 

Also  in  this  chapter,  part  D  describes  and  derives  the  Cramer-Rao 
bound  for  the  estimated  fo.  Part  E  extends  the  method  to  include  an 
arbitrary  signal  noise  power  function.  Part  F  discusses  a  method  for 
obtaining  initial  estimates  of  the  and  parameters  which  can  then  be 
used  in  the  estimation  procedure  of  part  C. 

B.  Preprocessing  of  the  Data 

In  the  absence  of  interchannel  coupling,  each  frequency  channel 
would  be  a  chi-square  distributed  random  variable  with  two  degrees  of 
freedom.  However,  since  each  filter  in  the  bank  of  filters  has  a  finite, 
but  unknown,  bandwidth  larger  than  the  separation  between  channels,  the 
received  record  is  then  expressed  as  the  original  power  spectrum  convolved 
with  some  "smearing  function".  Consequently,  the  degrees  of  freedom  of 
each  channel  Is  increased. 

The  estimation  method  requires  that  the  distribution  of  the  data  be 
of  known  form.  Furthermore,  computational  efficiency  is  gained  if  the 
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channels  are  uncorrelated  among  themselves.  The  solution  is  to  remove  the 
effect  of  coupling  by  filtering  the  record  in  the  frequency  domain  by  a 
digital  filter: 


X(fk) 


S(V,)  h,(m) 


k-1,2 . N 


X(fk)  is  a  nearly-Gausslan,  zero  mean,  uncorrelated  (and  therefore  independ¬ 
ent)  process.  The  filter  hj(m)  is  a  nonrecursive  zero-phase  (s3nnmetric) 
digital  filter,  designed  by  the  McClellan-Parks  algorithm  [1],  of  length 
2M+1  =  41  with  an  attenuation  of  40  dB  in  the  stopband.  The  problem 
encountered  at  the  edges  of  the  record  is  dealt  with  by  extending  the 
record  at  each  end  by  its  mirror  image 


^  “X  ®  5  ('N-n) 

for  n  =  1,...,M  and  N  =  128. 


This  processed  record  can  then  be  expressed  as 


X(fk)  =  N(fk)  +  H^(fk)  T(fk) 


where  NCf^)  and  TCf^)  are  Independent  Gaussian  variables,  with  zero  mean, 

but  a  variance  which  depends  on  the  signal  strength  S„  and  T„.  The 

N  N 

variance  of  ^Cf^)  is 
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C.  Maximum  Likelihood  Estimation 

Parameter  estimation  via  the  maximum  likelihood  technique  can  be 
expressed  intuitively  as  follows.  We  have  a  model,  with  parameters  that 
can  be  chosen.  Given  an  observation  X,  what  choice  of  parameters  make 
the  given  observation  most  likely  to  have  happened?  In  the  context  of  the 
above  model,  the  parameter  vector  is 


e  =■ 


] 


and  the  observation  vector  is  the  preprocessed  record 


X  =  [X(f^),  X(f^),  ....  X(fjj)]  *  [x^,X2,...,Xjj]. 

The  (ML)  estimation  procedure  Involves  finding  the  probability 
density  of  X  as  a  function  of  the  parameter  vector  0,  p(X;0)  which  Is  then 
maximized  with  respect  to  0.  In  practice,  the  likelihood  function  1(0) 
is  maximized,  where 


1(0)  -  log  p(X;0). 


After  preprocessing,  the  channels  are  Independent  zero-mean  Gaussian 
variables.  The  probability  density  can  be  expressed  simply  as 


p  (X;  -  n 
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exp 
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where  o  is  obtained  from  (A). 
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The  likelihood  function  Is 
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The  first  term  can  be  ignored  since  it  does  not  depend  on  the  parameters 
The  ML  estimate  is  obtained  from 


maximize  1(6) 
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E 
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log  a 


We  can  try  to  locate  the  maximum  by  taking  the  partial  derivatives 
2  2 

with  respect  to  f^,  and  eauate  to  zero.  However,  this  yields  a 

complex  set  of  coupled  nonlinear  equations,  and  the  solution  of  these 
equations  does  not  guarantee  a  global  maximum. 

Instead,  the  solution  is  found  numerically.  This  is  a  three- 
dimensional  maximization,  and  some  gradient  directed  algorithm  would  seem 
to  be  called  for.  However,  1(9)  is  itself  a  random  process,  making  its 
gradient  too  unreliable  for  such  a  sensitive  algorithm.  Instead,  the 
following  iterative  approach  is  used: 

2  2 

1.  Obtain  initial  estimates  for  and  .  These  can  be  rough 

guesses,  or  they  can  be  obtained  by  the  method  in  part  F. 

2  2 

2.  Keeping  Oj^  and  constant,  maximize  1  by  a  one-dimensional 

Fibonacci  search  on  f  . 


2  2 
3.  Keeping  O  and  f  constant,  do  Fibonacci  search  on  . 

o  2  T 

Keeping  and  constant,  search  for  . 

5.  Go  back  to  step  2  and  repeat  until  convergence  is  achieved. 

This  method  has  been  shown  to  be  reliable  with  the  given  data,  and 

converges  in  about  three  or  four  iterations,  when  initial  estimates  of 
2  2 

and  are  provided  by  the  method  of  part  F. 

N  T 

D.  Cramer-Rao  Bound 

Since  the  estimated  f  is  a  random  variable,  it  would  be  useful  to 

o 

determine  its  mean  and  variance,  and  hence  determine  the  accuracy  of  the 
estimator.  For  most  nonlinear  estimation  problems,  this  is  a  difficult 
task.  However,  the  theory  does  provide  a  quantity  known  as  the  Cramer- 
Rao  bound,  which  is  a  lower  bound  on  the  variance  of  the  estimate.  This 
bound  depends  only  on  the  model  and  the  statistics  of  the  data,  and  is 
independent  of  the  actual  estimation  technique  used.  While  it  is  not  as 
useful  as  an  upper  bound,  it  is  an  inherent  property  of  the  data  Itself, 
and  can  tell  us  if  our  expectations  concerning  estimation  accuracy  are 
reasonable . 

There  are  two  formulas  for  calculating  this  bound  [2].  In  this  case, 
the  easier  one  is 
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We  are  Interested  only  in  the  parameter  f^.  Suppose  and  are 

constants.  The  likelihood  function  is 
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where 


and  define 
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Taking  the  expectation  over  X,  ^  *^x  ^^k^  results  in 
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and  the  bound  is 
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2  2 

The  bound  depends  on  the  parameters  f  ,  o„  ,  as  well  as  on  the 

O  N  1 

filter  order  n.  Sufficient  accuracy  can  be  obtained  by  using  the  esti¬ 
mated  values  of  the  parameters  instead  of  the  true  values.  Generally,  the 

bound  is  lower  for  higher  values  of  filter  order  n.  It  also  decreases  as 
2 

Ojj"  (background  noise  power)  decreases. 

Although  there  is  no  guarantee  that  an  estimator  meets  the  Cramer- 
Rao  bound,  some  good  properties  of  the  ML  estimator  recommend  its  use. 

The  ML  estimator  is  consistent,  l.e.,  as  N  (the  number  of  observations)  in¬ 
creases.  the  estimate  converges  to  the  true  value.  The  ML  estimator  is 
also  asymptotically  efficient,  l.e.,  as  N  Increases,  the  variance  tends 
toward  the  Cramer-Pao  bound.  Thus,  the  performance  of  a  ML  estimator  can 
be  Improved  by  sampling  the  frequency  spectrum  more  finely. 


E.  Extension  to  Arbitrary  S^(f) 


In  the  previous  sections,  we  assume  that  the  variance  of  T(f)  is 
constant  over  frequency.  This  assumption  can  be  relaxed  if  information 
about  the  transmitted  terrestrial  radio  noise  is  available.  Then  equa¬ 
tion  (3)  can  be  modified  as 


X(fk)  =  N(f^)  + 


'  ^o  <^k> 


where  the  variance  of  T  dos  not  depei.d  on  frequency.  R(f,  )  is  a  function 

o  K 

which  expresses  the  relative  terrestrial  noise  power  in  the  frequency  band. 
It  can  be  multiplied  by  any  scale  factor  as  convenient.  The  variance  of 
X  becomes 


<'k> 
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which  can  then  be  used  in  the  maximum  likelihood  formulation  of  equation 
(7). 


The  Cramer-Rao  bound  becomes 
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Thus  If  the  signal  power  being  transmitted  near  the  subsatellite 
point  is  known,  it  can  be  Included  in  the  above  technique.  In  the  absence 
of  such  Information,  the  assumption  “  constant  must  be  made. 

2  2 

F.  Initial  Estimates  of  oN  and  C5T 


Estimates  of  and  are  obtained  by  another  maximum  likelihood 
technique  used  in  supervised  learning  [3].  Suppose  the  ionosphere  be¬ 
haves  as  an  ideal  filter.  Then  there  are  two  cases 


Wl;  Case  1.  only  present  =*  =  °1^ 

W2:  Case  2.  Sjj(f)  and  S^(f)  present  +*^1^  «  2^ 

The  probability  density  of  Xj^  can  be  described  as  the  following 
mixture  density 

^  2 
p(Xj^;0)  -  p(Xj^  ^ 

where  P(Wj^)  is  the  proportion  of  channels  for  which  Case  1  applies,  and 

P(W  )  the  proportion  for  which  Case  2  applies,  with  P(W  )  +  P(W  )  -  1, 
l  2  2^^ 

and  the  parameter  vector  is  £  is  (  ,02).  The  log  likelihood  func¬ 

tion  is 

N 

1  (0)  »  E  log  p(Xj^;£)  (13) 
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ttie  derivative  with  respect  to  is 
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Noting  that  \  :  =  p  (X^  Wj  ,  )  P(W^) 
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we  can  obtain  the  useful  form 
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Now  since  Is  Gaussian  and  zero  mean, 
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and  then 
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The  above  cannot  be  solved  in  closed  form,  but  an  Iterative  approach 


can  be  used: 


where 
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The  P  (W^)  are  also  unknown.  They  are  also  found  iteratively  by 

P(V.j)  "  FWjlSc:  <^1^)  (18) 

K*i 

Equations  (16-18)  provide  the  method  used  to  obtain  Initial  estimates 
2  2 

of  the  ,  and  therefore  .  Typically  this  algorithm  is  run  for  about 
ten  iterations,  and  the  results  then  used  In  the  estimation  method  of  Part 
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OPAQUE  -  A  MEASUREMENT  PROGRAM  ON  OPTICAL  ATMOSPHERIC  QUANTITIES 

IN  EUROPE 

PRuGRAM  PURPOSE  AND  OBJECTIV i 

The  ultimate  objective  for  the  research  program  Is  to  develop  a  data 
base  of  atmospheric  optical  and  Infrared  (IR)  parameters  which  affect  military 
systems.  Until  the  present  time,  efforts  to  measure  atmospheric  optical 
parameters  have  been  on  a  purely  national  basis.  The  time  has  now  come  to 
bring  those  efforts  together  and  to  consolidate  the  results,  thus  provlefliig 
a  more  universally  acceptable  data  base.  Possible  applications  for  the 
resulting  data,  would  be  In  the  optimization  of  optical  or  IR  vision  systems 
by  taking  due  account  of  the  influence  of  the  atmosphere  on  the  system 
capability.  For  example,  one  can  immediately  envisage  the  data  being  used 
in  the  optimization  of  Night  Recon,  Weapon  Guidance  and  Laser  Systems.  The 
results  from  the  measurement  program  would  also  be  useful  for  specifying  the 
atmospheric  optical  and  IR  (Infrared)  properties  for  operational  (system 
deployment)  purposes. 

It  will  be  possible  to  derive  from  the  measured  data,  distribution  of 
the  frequency  of  occurrence  for  various  parameters  such  as  Illumination 
levels,  transmission,  contrast  reduction,  as  a  function  of  geographical 
location,  time,  wave-length  region,  altitude,  etc.  Attempts  are  also  being 
made  to  derive  from  the  statistics,  material  relationships  between  the 
atmospheric  optical  parameters  themselves  as  well  as  with  the  general 
meteorological  condtlons,  which  would  be  useful  for  prediction  purposes. 

Electromagnetic  radiation  from  the  ultraviolet  through  the  visible 
and  Infrared  spectrum,  Is  affected  In  Its  propagation  through  the  earth's 
atmosphere  by  absorption  and  scattering  from  air  molecules  and  aerosol  or 
haze  particles  and  by  the  refractive  effects  of  atmospheric  density  changes. 
Therefore,  any  optical  system  whose  function  It  Is  to  look  through  the 
atmosphere  Is  affected  in  its  "seeing"  capabilities  by  these  atmospheric 
effects.  Optical  systems  are  being  used  in  such  military  functions  as 
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reconnaissance,  target  detection,  recognition,  and  acquisition;  also,  in 
ranging  and  target  designation,  in  weapons  guidance  and  in  optical  communica¬ 
tions,  Just  to  mention  some  examples.  Optical  systems  use  a  wide  range  of 
wavelengths  from  the  pho topic  response  of  the  human  eye  through  narrow  and 
broadband  systems  in  the  visible  and  infrared  (IR)  and  quasimonochromatic 
laser  systems. 

The  spectral  band  width  of  such  optical  systems  generally  has  to  be 
a  compromise  between  the  desire  to  keep  it  as  narrow  as  possible  for  reject¬ 
ion  of  undesirable  background,  and  the  need  to  make  it  wider  to  be  able  to 
receive  sufficient  radiation  energy  on  the  detector.  Active  laser  systems, 
on  the  other  hand,  have  created  now  the  possibility  for  extremely  narrow 
path  band  systems. 

Atmospheric  optical  properties  affect  systems  in  different  ways. 

In  all  cases  atmospheric  extinction  reduces  the  optical  signal,  originating 
from  the  object  scene,  along  its  path  through  the  atmosphere  to  the  receiver. 
This  atmospheric  extinction  is  due  to  the  loss  of  light  energy  in  the  direct 
beam  caused  by  molecular  and  aerosol  absorption  and  also  due  to  light 
scattered  out  of  the  direction  of  the  beam  by  molecules  and  aerosols. 

Some  systems,  such  as  reconnaissance  systems  or  contrast  seekers  for 
weapons  guidance,  depend  on  the  detection  of  brightness  differences,  or  con¬ 
trast,  between  object  scene  elements;  for  Instance,  a  target  and  its 
surrounding  background.  The  inherent  contrast  between  two  different  elements 
in  an  object  scene  is  reduced  by  the  atmosphere  because  of  light  scattered 
along  the  path  into  the  direction  of  the  receiver.  This  scattered  light  is 
called  air  light  or  path  radiance.  If  the  Intensity  of  path  radiance 
becomes  of  the  same  magnitude  as  the  radiance  from  the  object  scene,  the  con¬ 
trast  begins  to  disappear  and  different  scene  elements  can  no  longer  be 
distinguished . 

Imaging  systems,  looking  over  very  long  paths,  require  very  high 
angular  resolution.  In  such  systems,  the  resolution  is  also  reduced  by 
atmospheric  turbulence  and  resulting  optical  scintillation  phenomena  along 
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the  path.  The  picture  becomes  blurred  and  sharp  edges  become  fuzzy. 

Many  optical  military  systems  use  active  illuminatiDn  systems,  in 
particular  lasers  (for  example,  target  designator  systems,  range  finders). 

In  these  kinds  of  systems,  one  is  concerned  with  the  possibility  that  the 
laser  beam  might  be  detected  from  off-axis  scattering  by  atmospheric  aerosols 
and  molecules. 

Many  studies  have  been  conducted  to  determine  the  magnitude  of  these 
various  atmospheric  effects  on  military  systems  and  operations.  For  example, 
in  the  employment  of  the  F4/MAVERICK/EQ  missile,  cloud  ceiling  is  a  determin¬ 
ing  parameter  for  aircraft  attack  profiles;  visibility  affects  the  lock-on 
range  for  the  target-background  contrast  sensor  of  the  Maverick  missile. 

In  another  area,  this  data  can  be  used  to  estimate  the  scattered 
light  flux  on  an  aircraft  borne  sensor  originating  from  the  scattering  of  a 
laser  beam  which  is  not  directly  hitting  the  aircraft. 

The  OPAQUE  program  is  operating  a  number  of  stations  covering  the 
areas  of  Scandinavia,  Central,  Western  and  Southern  Europe  including  some 
island  or  shipboard  stations. 

The  choice  for  the  location  of  these  sites  was  based  first  on  military 
tactical  considerations;  secondly  on  their  representativeness  with  respect 
to  the  geographical  environment  and  atmospheric  meterological  condtions; 
and  thirdly  on  availability  of  logistical  support. 

All  data  is  being  reduced  to  a  commonly  agreed  upon  data  format 
compatible  with  computerized  data  analysis.  All  the  data  is  to  be  reduced 
to  an  agreed  format  on  a  six  monthly  basis,  and  then  collected  in  a  central 
data  bank  from  which  the  participating  nations  will  be  able  to  obtain  con¬ 
solidated  data  tapes.  Distribution  of  the  consolidated  data  will  be  at  an 
agreed  classification  level.  Each  member  nation,  or  group  of  member  nations 
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may  then  perform  their  own  anaJysls.  Published  reports  will  then  be  sub¬ 
mitted  to  a  central  secretary,  to  provide  a  common  series  of  OPAQUE  reports 
for  distribution  to  the  participating  nations  and  to  others  as  may  be  agreed 
iointly.  It  Is  also  desirable  that  a  six  monthly  review  be  made  of  the 
data  analysis  with  the  aim  of  avoiding  unnecessary  duplication. 


PHYSCIAL  BACKGROUND 

Every  optical  of  IR  system  which  is  looking  through  the  atmosphere, 
is  affected  by  atmospheric  light  scattering,  absorption,  turbulence,  atmosph¬ 
eric  refraction,  and  nonlinear  effects.  The  Important  parameters  for  any 
given  line  of  sight  for  the  atmospheric  transmission,  the  atmospheric  path 
radiance,  and  the  turbulence  factor. 

These  quantities  are  functions  of  the  illumination  conditions,  the 
atmospheric  properties  themselves,  and  also  the  average  terrain-background 
albedo.  These  factors  affect  the  amount  of  radiation  energy  transmitted 
from  a  source  along  an  atmospheric  path  to  a  receiver,  the  reduction  of  Image 
constrast,  and  the  distribution  (spatial  and  temporal)  of  energy  in  the 
transmitted  light  beam. 

lor  example,  under  homogeneous  conditions  the  apparent  luminance  of 
target  and  background  at  range  R  can  be  written  as 
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defining  the  contrast  at  zero  range  as 
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and  the  apparent  contrast  at  range  R  as 
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it  can  be  shown  that 
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Alternatively,  using  the  diffuse  reflectance  of  the  target  and 
background  A  and  A',  Eq.  (1)  can  be  written  as 

C 

C  =  _ 5 _ _ _ 

1  +  - - (e  R  -  1) 

(A  +  A’)E 


The  following  definitions  are  relevant  to  the  previous  equations 


(1)  C 

o 

(2)  L.L' 

(3)  T 

(4)  o 

(5)  R 

(6)  E 
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(7)  A, A’ 


=  Contrast  at  zero  range 

-2 

=  Luminance  of  target  and  background  [Cd  M  ] 

“  Transmission  =  e 
“  Attenuation  coefficient  [M 
“  Range  (m) 

=  Illuminance  indident  on  a  vertical  plane  [Lux] 
=  Diffuse  reflectance  of  target  and  background 


(3) 

-  Horizon  luminance  [Cd  M  ] 

(9) 
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o 

=  1  sr 

(10) 

Lp 

“2 

*  Path  Luminance  (Cd  M  )  per  range  R 

From  a  study  of  these  parameters  and  Eqs.  (1)  and  (2),  It  Is  possible 
to  specify  a  minimum  number  of  parameters  which  must  be  measured  in  any 
program  before  worthwhile  results  can  be  obtained.  Section  3  gives  what  are 
considered  to  be  the  minimum  set  which  It  Is  worthwhile  to  measure,  Including 
IR  transmission. 

is  the  off-axis  angle  (0°  would  be  the  direct  beam).  Such  signals 
could  be  detected,  and  in  response  an  aircraft  could  take  appropriate  counter¬ 
actions  . 


The  effect  of  illumination  levels  on  the  ability  of  an  observer  to 
detect  targets  of  various  sizes  has  gained  interest  and  lately  with  the 
development  of  image  intensifier  systems  and  Low  Light  Level  Television 
(LLLTV) .  These  systems  have  extended  passive  visual  seeing  from  twilight  down 
to  nighttime  illumination  levels. 

The  aerosol  scattering  and  absorption  properties  are  described  in 
the  Mie  theory.  This  theory  applies  to  spherical  particles  and  a  few  other 
simple  particle  configurations;  however,  comparisons  with  experimental 
results  show  good  enough  agreement  to  justify  for  most  practical  purposes 
the  application  of  the  Mie  theory  to  natural  aerosols.  The  scattering  and 
absorption  coefficients  for  aerosol  particles  between  sizes  r^^  and  and 
distribution  N(r)  are 


53 


and 


k 


Aa 


N(r)  dr 


Q  and  Q  are  functions  of  the  ratio  of  particle  size  to  wavelength,  and  the 
aerosol  particle  refractive  index  m;  they  are  derived  from  the  Mie  theory 
and  are  called  efficiency  factors  for  scattering  and  absorption,  respectively. 
The  refractive  index  of  aerosol  particles  is  a  real  number  of  completely 
transparent  substances  (no  absorption);  for  absorbing  substances  it  is  a 
complex  number. 


m  =  n  -  in ' 

where  n  is  the  real  part  of  the  refractive  index  and  n'  the  imaginary  part. 
The  quantity  Arm' /A  is  the  absorption  coefficient  of  the  aerosol  substance; 
its  dimension  is  (cm) 

These  equations  define  all  the  aerosol  absorption  and  scattering 
quantities  which  are  needed  to  calculate  the  transmission  of  a  light  beam 
through  the  atmosphere. 

The  Raleigh  and  Mie  theory  also  give  relationships  for  the  intensity 
of  radiation  scattered  out  of  a  light  beam  into  different  directions.  For 
molecular  Rayleigh  scattering,  the  angular  scattering  intensities  are  given 
by 


\  <»> 


(1  +  cos 


cp) 
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For  spherical  particles  with  sizes  of  the  same  order  of  magnitude 
or  larger  than  the  wavelength  of  the  incident  light,  the  Mle  theory  gives 
for  the  angular  scattering  Intensltlties 


■x  <®>  ■  \  o 


1 
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In  order  to  predict  the  atmospheric  optlcal/IR  effects  on  the 
performance  of  a  given  optical  system,  two  things  art  needed:  (a)  algorithms 
to  calculate  the  specific  quantity  required  for  the  physical  condtlons  under 
which  a  system  functions  (for  example,  the  transmittance  for  a  specified 
spectral  band  and  along  a  given  atmospheric  slant  path) ,  and  (b)  a 
description  of  the  atmospheric  properties  which  are  required  as  input  to 
(a)  (for  example  the  aerosol  and  molecular  extinction  coefficients  as  a 
function  of  wavelength,  and  the  vertical  profile  of  these  parameters).  The 
development  of  the  algorithms  is  primarily  a  theoretical  effort;  the  atmosph¬ 
eric  properties  must  be  derived  from  field  measurements. 


ATMOSPHERIC  DATA  BASIS 

The  basic  atmospheric  quantities  which  are  needed  to  calculate 
atmospheric  optical  effects  are:  the  distribution  of  atmospheric  molecules 
and  aerosols  and  their  absorption  and  scattering  properties,  and  the  turbu¬ 
lence  structure  constant  C  1  For  most  of  these  parameters,  models  describing 

n 

their  spatial  distribution  and  in  some  cases  temporal  variations  have  been 
developed . 

Models  for  the  distribution  of  atmospheric  molecular  density,  tempera¬ 
ture  and  water  vapor  concentration  are  defined  in  the  "Standard  Atmospheres" 
and  are  based  on  long  series  of  measurements.  The  data  base  for  other  minor 
gaseous  components  are  less  extensive.  Laboratory  and  field  measurements 


55 


exist  of  optical  properties  of  the  atmospheric  molecular  compnents,  and  have 
been  compiled  for  optical  calculations.  Data  on  over  100,000  absorption 
lines  from  different  molecular  species  have  been  compiled  for  wavelengths 
from  the  visible  into  the  far  infrared. 

Whereas  most  of  the  molecular  absorption  is  caused  by  two  or  three 
mo'ccular  species  only  (water  vapor,  carbon  dioxide  and  ozone),  the  extinct¬ 
ion  due  to  atmospheric  particulates  or  aerosol  particles  is  much  more  complex. 
The  concentration,  size  distribution,  and  composition  of  aerosols  is 
extremely  variable,  spatially  and  in  time;  the  development  of  representative 
models  for  the  aerosol  properties,  therefore,  is  in  a  more  preliminary  state. 
Aerosol  models  have  been  developed  for  continental  environments,  such  as 
rural  or  urban  areas,  and  for  maritime  regions.  These  models  describe  also 
the  variation  of  aerosol  properties  with  altitude  from  the  surface  up  to 
100  km  altitude.  The  vicissitude  even  within  these  models,  however,  makes 
their  applicability  to  a  given  real  world  situation  very  difficult. 

With  the  exception  of  those  parameters,  which  are  being  collected 
routinely  as  part  of  the  standard  meterological  observations  (that  is, 
temperature,  humidity,  surface  visibility),  the  present  data  base  for 
atmospheric  optical  modeling  has  no  statistical  significance.  It  describes 
only  "general  average"  conditions. 

ALGORITHMS  FOR  OPTICAL  CALCULATIONS 

Model  calculations  can  be  only  as  good  as  the  input  data  for  them. 

It  can  be  stated  at  the  outset,  that  in  general  the  theoretical  concepts  and 
the  algorithms  for  computing  optical  propagation  characteristics  are 
superior  to  the  accuracy  of  the  input  data.  This  is  certainly  correct  for 
transmission,  emission,  and  single  scattering  computations,  and  even  for 
most  multiple  scattering  computations. 

Algorithms  for  predicting  scattering  intensities  are  needed  primarily 
for  contrast  reduction  and  sky  radiance  calculations.  Path  radiance  P  is 
composed  of  light  which  has  been  scattered  several  times  along  its  path 
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through  the  atmosphere.  Depending  on  how  optically  thick  the  atmosphere  Is, 
’ight  which  has  been  scattered  up  to  ten  times  or  more  may  contribute 
significantly  to  the  radiation  intensities  of  concern.  No  rigorous  analyti¬ 
cal  solution  has  been  presented  as  of  this  date  for  the  radiation  transfer 
problem  involving  multiple  scattering  in  a  real,  aerosol-containing  atmos- 
pliere . 

NEKD  FOR  EXTENDED  ATMOSPHERIC  OPTICAL  DATA  BASE 

The  military  community  is  confronted  with  two  different  types  of 
problems  in  predicting  the  atmospheric  optical  effects  on  the  performance 
of  n  given  optical  system.  The  first  one  is  associated  with  the  question: 
"Wliat  is  the  probability  that  a  given  system  will  be  able  to  function 
successfully  in  a  certain  environment?"  This  question  is  being  asked  by 
the  systems  designer  to  derive  the  appropriate  systems  specifications.  It 
is  also  being  asked  by  the  systems  analyst  to  develop  the  proper  operational 
procedures  and  deployment  plans  for  a  system. 

The  second  type  of  requirement  is  for  a  specific  forecast  for  the 
atmospheric  optical  environment  at  the  time  and  place  of  deployment  of  the 
system. 


The  first  problem  is  one  of  developing  an  atmospheric  optical 
climatology,  and  the  second  one  requires  a  capability  to  forecase  atmospheric 
optical  conditions  from  the  available  standard  meteorological  observations 
and  the  general  weather  forecast.  At  present,  however,  one  has  in  many 
cases  only  a  qualitative  understanding  of  the  physical  relationships  between 
the  meteorological  factors  and  the  environment  on  one  side,  and  the 
atmospheric  optical  properties  on  the  other.  A  quantitiative  definition  of 
these  relationships  requires  also  a  reliable  and  comprehensive  data  base  as 
a  source  for  empirical  correlations. 

While  such  data  are  useful  for  system  design  planning,  they  must  be 
used  cautiously  in  operational  planning.  The  data  will  be  misused  if 
employed  to  estimate  the  chances  that  a  single  or  small  number  of  missions 
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will  find  suitable  conditions  of  cloud  cover  and  surface  visibility.  If  the 
data  are  used  in  relation  to  large  number  of  missions,  then  more  conflcence 
can  be  attached  to  them. 

These  data  can  be  of  use  to  the  operational  weather  forecaster, 
however,  by  providing  him  useful  relationships  for  shaping  his  forecasts. 
Since  the  type  of  air  mass  is  so  important  in  influencing  visibility  and 
cloud  cover,  a  refinement  of  the  data  by  air  mass  type  of  synoptic  situation 
would  be  even  more  useful  to  the  forecaster. 
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The  atmospheric  environment  Is  of  a  quite  decisive  Influence  on  many 
electro-optical  systems.  In  moat  cases,  the  atmospheric  medium  constitutes 
a  limiting  factor  for  the  propagation  of  electro-magnetic  signals.  In  some 
other  cases,  there  are  the  optical  properties  of  this  medium  which  render 
the  systems  Infeasible.  Thus  to  support  electro-optics  system  planners, 
designers  and  users  with  probability  of  occurrence  Information  on  earth 
surface  propagation  conditions.  It  Is  desirable  to  have  statistics  of 
specific  propagation  parameters  available. 

The  natural  aerosol  Is  a  basic  constituent  of  the  atmosphere,  having 
considerable  impact  generally  on  optical  propagation  properties.  Its  basic 
parameters  are  composition,  concentration,  and  distribution. 

Though  considerable  effort  has  been  addressed  in  the  past  by  many 
research  workers  to  measuring  and  understanding  the  characteristics  of 
natural  aerosols,  the  state-of-the-art  must  be  considered  still  yet  as 
quite  preliminary.  In  fact  there  are  neither  comprehensive  statistics 
available  containing  sufficient  Information  for  purposes  mentioned  above, 
nor  are  there  unique  models  which  can  be  used  for  operational  purposes  free 
of  ambiguities. 

The  reason  for  this  lies  in  the  complexity  of  aerosol  composition,  the 
vast  mechanisms  of  its  formation,  Its  complex  dynamics,  and  In  Its  complex 
way  of  ageing.  One  other  problem  in  the  past  has  been  the  state-of-the-art 
of  data  collecting  devices  which  did  not  easily  permit  measuring  aerosol  data 
over  larger  size  intervals  with  satisfactory  resolution. 

The  situation  Is  improving  now  with  the  advent  of  new  aerosol  particle 
counters , 

Typically,  aerosol  size  distribution  measurements  generate  a  huge  stream 
of  original  data.  These  must  be  channeled  down,  mainly  for  two  reasons. 
First,  because  of  economics  all  redundancy  in  the  data  should  be  removed  to 
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the  largest  extent  possible.  Secondly,  there  is  virtually  no  way  to 
sensibly  work  with,  analyze,  or  use  in  later  processing  applications  this 
bulk  of  original  data. 

Generally,  in  a  situation  like  this  the  benefit  of  a  model  consists  in 
providing  some  kind  of  mathematical  expression  through  which  the  data  could 
be  greatly  compressed.  However,  there  is  no  such  single  expression  capable 
of  describing  size  distributions.  Junge’s  formula  is  referred  to  very 
widely  in  the  literature.  It  is  an  empirical  model  in  the  sense  that  observa¬ 
tions  show  that  it  tends  to  be  followed,  but  not  that  it  must  be  followed. 

Yet,  for  the  applications  metnloned  above,  Junge's  formula  does  not  seem  to 
be  adequate  enough.  Other  distribution  functions  are  the  standard  and  the 
modified  gamma  distribution  or  the  normal  and  the  lognormal  distributions. 
Gamma,  normal  and  lognormal  distributions  have  proven  to  be  quite  versatile. 
Generally  they  are  used  for  typical  background  aerosols.  Sums  of  normal  and 
lognormal  distributions  are  used  to  fit  distributions  originating  from  a 
single  source  through  a  single  mechanism.  These  few  examples  were  cited  to 
show  that  very  often  the  choice  of  a  specific  mathematical  expression  seems 
to  be  justified  only  by  the  purpose  of  use.  A  highly  accurate  method  for 
fitting  aerosol  data  is  available  through  the  use  of  spline  functions. 

However,  this  approach  normally  yields  a  large  number  of  coefficients  which 
do  not  correspond  to  sphysically  meaningful  parameters. 


2.  AEROSOL  SIZE  DISTRIBUTION  MEASUREMENTS 
2 . 1  Size  Distribution  Function 

Assume  an  aerosol  is  contained  in  a  unit  probing  volume.  The  basic 
method  for  determining  its  size  distribution  function  is  to  count  the  total 
number,  N(D)  of  aerosol  particles  (per  unit  volume)  whose  diameter  values 
fall  below  a  certain  value,  D.  The  function  N(D)  is  a  cumulative  function 
with  values  between  zero  and  N,  the  total  number  of  aerosol  particles  contain¬ 
ed  within  the  unit  volume.  Generally,  the  number  of  aerosol  particles  in  a 
given  volume  is  very  large,  and  their  sizes  may  be  assumed  to  vary  continuous¬ 
ly.  Thus,  N(D)  approaches  a  continuous,  smooth  function  as  N  Increases.  It 
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is  then  possible  to  derive  the  size  distribution  function,  n(D),  by 


„  (W  . 


(2.1) 


or  in  integral  representation, 


U 

J  ” 


(D')  dD’ 


(2.2) 


In  practice  the  determination  of  N(D)  leads  to  serious  problems  since 
it  is,  for  experimental  reasons,  possible  only  to  measure  particles  with 
sizes  above  a  certain,  device-dependent,  lower  limit,  say  D^.  Thus  N(D)  is 
determined  up  to  a  constant  N(Dj^).  Since  n(D)  is  a  derivative  function,  this 
ambiguity  is  removed. 

The  aerosol  counters  employed  in  this  work  provide  data  on  the  number 
of  particles  (per  unit  volume)  for  some  specified  diameter  intervals 

£  D  with  i  =  1,2,...,M.  If  denotes  the  number  of  particles 

for  the  i-th  channel,  then  because  of  (2.2)  one  has 


n(D’)  d  D' 


(2.3) 


Accordingly,  may  be  used  directly  to  evaluate  an  estimate  for  the  size 
distribution  function,  n,  at  the  point  that  is 


n(D^)  ~  n^:  = 


!  1  “  1,2,... ,M. 


(2.4) 


°1+1  -  °1 


61 


3. 


DATA  FIT  OBJECTIVES 


If  data  are  going  to  be  fitted  by  some  algorithm,  two  basically 
different  cases  are  to  be  distinguished.  In  the  one  case,  the  data  result 
from  a  process  the  physical  nature  of  which  is  understood  in  the  sense  that 
a  physical  model  or  theory  is  available.  In  this  case,  data  fitting  permits 
the  determination  of  some  free  model  parameters  through  adjusting  the  model 
predictions  to  the  set  of  measured  data.  In  the  other  case,  there  is  no 
physical  model  at  hand  to  predict  the  experimental  findings,  and  data  fits 
in  this  case  are  primarily  another  way  of  data  presentation.  If  the  fit 
algorithm  is  chosen  carefully  and  appropriately  adapted  to  the  problem,  a 
few  Important  objectives  may  be  met  simultaneously: 

A.  Data  representation: 

The  Insertion  of  a  smooth  fit  curve  into  the  plot  of  data  points  is 
a  useful  way  to  assist  the  understanding  of  the  general  trends 
underlying  the  data  set  by  inspection. 

B.  Error  smoothing: 

The  smooth  fit  curve  is  likely  to  smooth  out  high  frequency 
contributions  within  the  set  of  data  due  to  random  error.  This  is 
especially  true  if  a  least-squares-f it  algorithm  is  applied. 

C.  Data  reduction: 

Through  data  fitting,  one  can  generally  reduce  the  large  amount  of 
data  to  a  few  functional  parameters  only.  Frome  these,  it  is 
possible  to  retrieve  the  original  data  (in  the  accuracy  desired), 
as  well  as  to  calculate  new  Interpolated  data. 

D.  Empirical  modelling: 

Generally,  it  is  much  easier  to  correlate  any  typical  condition  of 
the  process  under  study  with  a  typical  pattern  of  the  reduced  data 
set,  than  with  the  original  samples. 
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The  Least  Squares  Method 

Measured  data  with  values  v^,  (  1  “  are  assumed  to  be  a 

function  of  a  single  variable,  t.  The  experimental  result  is  a  list  of  M 
data  pairs  (t^,v^).  The  measuring  points,  t^,  need  not  be  spaced  regularly. 
In  the  least-squares  method,  a  function  F  from  an  appropriate  collection  F 
is  determined  such  that 


S  (F(t  )  -  V  )  -  min. 

i=l 


(A.l) 


For  ease  of  computation,  the  function  space  F  is  chosen  as  to  be  formed  by 
all  linear  combinations  of  some  basis  functions,  fj^.  That  is,  any  F  eF  can 


be  written  as 


F(t)  -  E  b  f  (t)  . 
k-0 


(4.2) 


The  set  of  f^^  must  be  linearly  independent.  Though  desirable  in  many  cases, 
the  fj^  need  not  necessarily  obey  some  orthogonality  condition.  In  most 
applications,  the  function  space  F  is  deliberately  restricted  by  requiring 


K  <  M. 


(4.3) 


The  free  amplitudes,  bj^,  in  (4.2)  are  determined  through  (4.1)  as  solutions 
of  the  Gaussian  normal  equations 


A  .  B  -  C, 


(4.4) 


where  A  is  a  positive  definite,  square  symmetric  matrix  with  elements 


E  f.  (t.)  f.  (t,).  (k,  1,  »  1....K).  (4.5) 

1=1  ^ 
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and  C  is  a  column  vector  with  elements 


M 

fj^(t^)  .  v^.  (k-l,...K).  (4.6) 

The  vector  B  of  elements  b,  is  the  solution. 

k 

The  problem  encountered  in  determining  solutions  of  Eq.  (4.4)  is  that 
the  matrix  A  usually  tends  to  be  ill-conditioned.  As  a  result,  numerical 
instabilities  arise.  The  problem  is  reduced  somewhat  if  an  appropriate 
function  system,  fj^,  is  used  for  expansion.  Furthermore,  it  is  advisable 
not  to  apply  the  genral  Gaussian  inversion  technique  to  deduce  the  solution, 
B,  rather  than  make  use  of  Cholesky's  method,  a  numerically  stable  matrix 
decomposition  technique  which  fully  exploits  the  symmetry  and  positive 
definiteness  of  matrix  A. 

4 . 2  Data  Transformation 

Size  distributions  of  natural  aerosols  tend  to  drop  off  in  magnitude 
over  several  decades  within  the  diameter  range  from  .1  to  30  m.  Therefore, 
the  fit  usually  is  applied  not  to  the  original  data  n^  but  rather  to  their 
logarithms 

v^  =  log  n^^,  (1=1,..., M).  (4.7) 


It  should  be  noted  that  his  is,  in  fact,  a  substantial  step  in  the 
data  analysis.  Indeed,  if  the  fitted  curves  are  going  to  be  employed  in 
any  optical  propagation  calculation,  taking  the  logarithms  is  a  necessary 
approach.  It  has  been  shown  that  the  results  of  such  calculations  depend 
to  a  very  large  extent  on  the  number  of  small  aerosol  particles  as  well  as 
on  that  of  larj^jr  particles,  although  the  latter  usually  are  much  less 
populated.  For  the  purpose  of  optical  calculations,  the  "goodness"  of  fit 
should  be  unifori.  •’ver  the  entire  spectrum  of  particle  diameters  considered 
regardless  of  their  magnitude  in  concentration.  This  is  quite  naturally 
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achieved  in  using  the  logarithmic  scale  since  then  the  least-squares  fit 
results  in  an  overall  minimized  relative  error  for  the  original  data. 


Typically  for  the  aerosol  counters  employed  in  this  type  of  work,  the 
widths  of  their  different  size  channels  generally  Increase  with  increasing 
diameter  values, 


Chebyshev  Polynomials 


The  success  of  the  fitting  procedure  depends  critically  on  the  set  of 
basis  functions,  f ,  to  be  use<d  for  expansion  of  the  data  points.  In  the 
present  work  it  is  suggested  that  for  fitting  natural  aerosol  size  distribu¬ 
tions,  the  use  of  Chebyshev  polynomials  is  highly  favorable.  In  this  case, 
formula  4.2  is  called  a  Chebyshev  expansion  and  is  of  the  form 


F(t)  =  b^  T^  +  b^  T^  +  ...  +  b^  T^. 


(4.10) 


This  saves  computing  time  considerably,  and  enables  an  easy  check  for 
convergence.  However,  it  is  then  more  difficult  to  interpret  the  results 
of  the  fit  in  a  straight  forward  and  general  way. 

Chebyshev  polynomials  are  well  known  from  many  textbooks  on  approximation 
theory  and  smoothing.  There,  one  usually  may  also  find  an  outline  of  their 
special  properties  which  makes  them  a  rather  unique  set  of  basis  functions 
for  use  in  numberlcal  fit  algorithms.  The  first  few  members  of  this  set 
are  given  by  (10) 


T^(t): 

-  1 

o 

Tj(t): 

T2(t): 

«  t 

-  -1 

+2t^ 

-  3 

T-(t): 

=  -3t 

+  4t 

3 

2  „  4 

T,(t): 

15(1): 

“  1 

=  5t 

-  8t  +  8t 

-20t^  +16t^ 
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The  first  polynomials  to  are  shown  graphically  in  Figure  2  for 
the  interval  -1  t  1.  From  this  figure  it  may  be  seen  that 

I  Tj^(t)  1  ;  |t)  f  1.  C4.12) 


The  definitions  (4.11)  indicate  that  the  polynomials  are  built  up 
from  either  even  or  odd  power  functions  in  t.  The  amplitudes  of  the  power 
functions  obviously  tend  to  be  very  large  in  magnitude,  and  it  is  because 
of  their  alternate  signs  that  the  stay  within  the  limits  (4.12).  The 
numerical  evaluation  of  series  like  this  may  generally  be  subject  to  serious 
computational  errors  for  computers  having  limited  word  lengths.  A  convenient 
way  out  of  this  problem  is  to  use  the  three-terms  recurrence  relation 


\+2 


2t 


\-fl 


T^(t) 


k=0,l,2,, 


(4,13) 


In  conjunction  with  the  basic  definition  (4.11)  for  T  and  T, . 

0  1 

Some  Examples 

Though  the  examples  chosen  are  rather  simple,  they  are  of  general 
Interest  in  the  study  of  aerosols.  They  are: 

-  The  Junge  size  distribution 

-  A  unimodal  lognormal  size  distribution 

-  The  haze  L  model 


The  Junge  size  distribution  is  of  the  form 


n  (D): 


where  the  prime  is  to  indicate  that,  since  n  is  a  function  of  diameter 
rather  than  of  radius,  the  concentration  C  differs  in  value  from  those  to 
be  found  in  conventional  tabulations. 
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On  a  log-log-scale,  n(D)  plots  as  a  straight  line.  In  this  case, 
the  Chebyshev  fit  results  In  a  linear  function 


F(t) 


1 

k-0 


b  +  b,  t  . 
o  1 


(4.15) 


The  interpretation  of  the  two  amplitudes  b  and  b,  is  most  conveniently 

O  1 

done  in  the  double  logarithmic  plot.  It  may  be  seen  that  b^  mainly 
determines  the  absolute  concentration  whll*  determines  the  slope. 


Junge  distribution  and  Chebyshev  approximations 


The  Junge  parameters,  C*  and  B  can  be  expressed  uniquely  in  terms  of 
the  amplitudes,  and  b^^. 


log  C  -  b  -  b, 
o  1 


log  (Djj  / 


B  +  1  = 


log  (Djj  / 


(4.16) 


From  this  It  Is  clear  that  the  amplitudes  b^^  are  dependent  on  the  low 

and  high  ends,  D  and  D  respectively,  of  the  measured  diameter  interval. 

L  h 

Thus,  when  comparing  Chebyshev  amplitudes  for  measurements  taken  with 

different  devices,  the  corresponding  values  for  D  and  D  must  be  taken  into 

L  H 

consideration. 


In  case  of  Junge  distributions,  we  can  avoid  this  unpleasant  device- 
dependence.  However,  the  problem  is  much  more  involved  for  more  general 
distributions  containing  many  more  terms.  We  hope  to  resolve  this  problem 
in  future  studies. 


The  next  example  is  that  of  a  unimodal 


n(D) 


2N' 


D. 


log 


6 


exp 


1^ 

2 


lognormal 


distribution 

,  D 

log  5 

_ m 

log  6  g 


Again,  the  prime  indicates  that  in  the  current  notation  the  parameter 
N'  is  different  in  value  from  that  in  the  literature  where  n  has  been  con¬ 
sidered  to  be  a  function  of  the  particle  radius.  The  example  illustrated 
is  an  average  si^e  distribution  for  maritime  aerosols,  and  taken  here  at 
the  experimental  meshpolnts  of  the  present  work.  The  median  number  diameter 
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accordingly  is  D  =  0.52uni;  the  logarithmic  standard  deviation  Is  t  =  4.8. 

m  g 

One  can  see  that  the  Chebyshev  fit  of  a  unlmodal  lognormal  distribution 
is  quite  excellent.  Moreover,  three  coefficients  are  sufficient  to  specify 
this  fit,  i.e.,  b^ i  b^i  and  as  one  would  expect  since  a  lognormal  distri¬ 
bution  is  represented  in  the  double  logarithmic  scale  by  a  parabola  of 
order  2. 


The  physical  parameters  N',D^  and  6^  can  be  expressed  uniquely  in 
terms  of  the  amplitudes  b^ ,  b^^,  and  b^. 

The  last  example  we  consider  is  that  of  Del rmend j Ian ' s  Haze  L  model 
distribution 


nCD)  “  a'  exp  (-b'  •  d'^  ) 


(4.18) 


with  the  following  parameter  values 


a 

b’ 


=  2.4878  10^  r  urn  ,  a  =  2 


=  10.6905 


r  ] 


6=5 


Unlike  the  two  other  examples,  the  distribution  does  not  in  the  log- 
log-scale  represent  a  simple  polynomial.  Thus,  the  Chebyshev  fit  in  this 
case  is  really  an  approximation. 

The  few  examples  presented  have  been  chosen  to  demonstrate  the  main 
advantages  and  shortcomings  of  fit  time  algorithms.  The  fit  was  shown  to 
be  versatile  in  the  sense  that  various  size  distributions  can  be  treated 
although  they  are  of  quite  different  nature.  All  of  them  are  basically 
related  to  the  aerosol  field.  Moreover,  unlike  in  cases  of  orthogonal 
polynomials,  the  fixed  basis  of  Chebyshev  interpretation  of  the  fitted 
data  in  terms  of  a  reduced  data  set,  namely  the  resulting  Chebyshev 
amplitudes.  A  remaining  difficulty  lies  in  the  fact  that  the  resulting 
amplitudes  are  dependent  on  the  range  of  particle  diameters  covered  in  the 

fit,  and  therefore  to  some  extent  are  device-dependent. 
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OPAQUE  DATA  BASE 


I.  DATA  BASE 

A.  Various  data  base  definitions  have  been  examined  and  implemented 

for  test  purposes.  Questions  in  this  area  revolve  around  two  points 

of  concern; 

1  .  Data  Structure  -  A  relatively  stable  logical  structure  has 

evolved  since  the  data  content  is  well  defined  and  understood, 

2.  Data  Retrieval  -  (Repeating  group  time  assignments,  element 

designation)  -  A  good  deal  of  experimentation  and  consideration 
has  been  given  to  the  method  of  storing  data  i.e.  the  representa¬ 
tion  of  a  logical  entry  (HR.,  DAY,  MO. 5.  This  will  play  a  key 
role  in  determining  data  base  size  and  retrieval  and  update 
capability.  Hour  as  the  logical  entry  will  allow  data  history 
storage,  while  day  as  logical  entry  requires  data  sets  to 
represent  hours  of  data. 

B.  Data  Base  Loading  -  Two  methods  of  initial  data  base  loading  have 
been  considered; 

1.  Natual  language  "load"  capability  -  natural  language  provides  a 
rudimentary  capability  to  mass  load  a  data  base  with  machine 
stored  data.  This  capability  has  been  examined  and  tested. 
Although  successful,  the  method  is  cumbersome  requiring  a  two- 
step  process. 

2.  PLI  Capability  -  PLI  provides  more  powerful  tools  for  data  base 
communication  i.e.  use  of  Host  Fortran  Processing)  and  appears 
to  be  the  method  which  will  be  chosen  for  data  base  loading. 
Preliminary  Fortran  coding  has  been  completed  and  PLI  statement 
requirements  examined,  but  no  PLI  tests  liave  been  performed, 

C.  Data  Base  Maintenance  -  The  following  areas  are  of  Interest; 

1.  Retrieval  -  plans  have  been  made  to  incorporate  PLI  statements 
into  existing  interactive  graphics  software.  No  testing  has 
occurred  but  preliminary  coding  has  been  completed. 
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2,  Update  -  a  procedure  has  been  defined  to  allow  Interactive 
editing  of  OPAQUE  data  resident  in  the  data  base. 

3.  DBA  Functions  -  a  procedure  to  provide  for  transaction  log, 
automatic  backup,  reco'-ery,  etc.  has  been  defined. 

II.  INTERACTIVE  PLOTTING  SOFTWARE 

A.  Interactive  plotting  software  currently  exists  and  will  provide 
single  plots  of  any  data  word.  No  provision  exists  to  test  quality 
value  or  data  out  of  range,  but  log  or  linear  plots  are  available, 
automatic  or  manual  7-Axls  scaling  Is  possible,  labels  are  provided, 
and  titles  may  be  entered  by  the  user. 

B.  A  second  version  uf  (A)  above  exists  which  Incorporates  those 
features  absent  from  (A)  but  no  testing  has  occurred.  This  system 
will  be  split  to  form  a  log  and  a  linear  program. 
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LOADN  77 


Step  1 


MET  DATA 
(3  MONTHS) 


LOWELL  FILE 
(1  MONTH) 


OP  FILE 


Step  2 


•PLOTTING 

PACKAGE 


Step  3 


UNCORRECTED 
OPAQUE  FILE  V 
“X  Step  4 


TOLERANCE  FOR 
.  80,  83  &  79 


QUALITY  UPDATE 


Step  5 


FINAL  TOLERANCES/  Step  6 
FOR  BAD  DATA  I 


LIST  BY 
INSTRUMENT 


Step  7 


FINAL 
COMPARE , 


Step  8 


PACKED 

DATA 
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array  specification  -  DIRFXT  ACCESS  FILE  -  DATA  BASE  SOFTWARE 


TYPE:  INTEGER 
DIMENSION;  140 
CONTENT  BY  WORD: 


ARRAY  WORD  if 

1 

2 

3 

4 

5 

6 

7 

8 
9 

10 

11 

12 

13 

14 

15 

16 

17 

18 

19 

20 
21 
22 

23 

24 

25 

26 

27 

28 

29 

30 

31 

32 

33 

34 

35 

36 


Station  No. 

Date 

Time 

Duration  of  Meas.  Cyc. 
Comment  Code  1 


OPAQUE  WORD  if 


Scattering  Value 
Filter  Value 
Humidity  Value 
MRI  Photopic,  Beg.  Val. 

"  "  ,  Fin.  Val. 

"  "  ,  Max.  Val. 

"  "  ,  Min.  Val. 

"  ,  No.  Meas. 

Eltro  Transmissoraeter ,  Beg,  Val, 
"  "  ,  Fin.  Val. 

"  "  ,  Max,  Val. 

"  "  ,  Min.  Val. 

"  "  ,  No.  Meas. 

Horizontal  Luxmeter,  Beg,  Val. 

"  "  ,  Fin.  Val. 

"  "  ,  Max.  Val. 

"  "  ,  Min.  Val. 

"  "  ,  No.  Meas. 

Vertical  Luxmeter,  North  Value 
"  "  ,  East  Val. 

"  "  ,  West  Val, 

"  "  ,  South  Val. 

Night  Path  Luminance,  Beg.  Val. 

"  "  "  ,  Fin,  Val. 

"  "  "  ,  Max.  Val. 

"  "  "  ,  Min.  Val. 

"  "  "  ,  No.  Meas. 


FORMAT 

12 

16 

14 

12 

13 

13 

13 

13 

13 
11 
11 
11 

14 
14 
14 
14 

13 

14 
14 
14 
14 

13 

14 
14 
14 
14 

13 

14 
14 
14 
14 
14 
14 
14 
14 
13 
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ARRAY  WORD  it 

DATA 

OPAQUE  WORD  If 

FORMAT 

37 

Variance  Path  Pet.  Meter,  Beg.  Val. 

35 

14 

38 

tl  |(  tl 

"  ,  Fin.  Val. 

36 

14 

39 

M  »l  tt 

’•  ,  Max.  Val. 

37 

14 

AO 

n  fi  It 

"  ,  Min.  Val. 

38 

14 

A1 

tl  tt  tl 

"  ,  No  Meas. 

39 

13 

42 

It  tt  tl 

"  ,  South  Val. 

40 

14 

A3 

tl  tl  It 

"  ,  West  Val. 

41 

14 

AA 

II  It  tt 

"  ,  North  Val. 

42 

14 

A5 

Filtered  Epplet  Data, 

-  .945 

43 

14 

A6 

It  tt  t* 

-  .4 

44 

14 

A7 

tl  It  It 

-  .87 

45 

14 

A8 

tt  II  tt 

“  1.06 

46 

14 

A9 

II  tl  tt 

=  .75 

47 

14 

50 

ft  It  rt 

=  .55 

48 

14 

51 

tl  It  tt 

Photopic 

49 

14 

52 

tt  tl  tt 

+  .3-3.51 

50 

14 

53 

It  It  tt 

Direct 

51 

14 

5A 

II  It  tt 

-  .3-3.52 

52 

14 

55 

Barnes  Trans.,  3-5 

53 

14 

56 

"  "  ,  8-12 

54 

14 

57 

"  "  ,  8-13 

55 

14 

58 

"  "  ,  Open  or  A 

56 

14 

59 

"  "  ,  3-5  Fin 

57 

14 

60 

Cloud  Cover  „ 

68 

11 

61 

Wind  Dir .  at  10 

69 

12 

62 

Wind  Speed  at  10 

70 

12 

63 

Wind  Dir.  at  2 

71 

12 

6A 

Wind  Speed  at  2 

72 

12 

65 

Pressure 

73 

13 

66 

Temperature 

74 

13 

67 

Dew  Point  Temperature 

75 

13 

68 

P«aln  Rare 

76 

13 

69 

General  Ground  State 

77 

11 

70 

Total  Rain  Past  Hour 

85 

13 

71 

Contel  and  Add.  Meas. 

Code 

64-1 

11 

72 

Aerosol  Ins.  Code 

58-1 

11 

73 

Aerosol  Date  Code 

58-2 

11 

7A 

Max.  Channel  Used 

58-4 

11 

75 

Aerosol  Rel .  Code 

58-3 

11 

76 

Add.  Extinction  Val. 

64-3,4,5 

13 

77 

Add.  Extinction  Val. 

Rel.  Code 

65-1 

11 

78 

C02  Concentration  Val 

• 

65-2,3,4 

13 

79 

"  "  Rel 

.  Code 

65-5 

11 

75 


ARRAY  WORD  # 

DATA 

OPAQUE  WORD  // 

FORMAT 

123 

Barnes 

Transmlssomecer 

Qual. 

.55 

84-1 

11 

124 

It 

tl 

It 

.56 

84-2 

11 

125 

II 

tl 

II 

,57 

84-3 

11 

126 

tl 

It 

ft 

.57 

84-4 

11 

127 

It 

It 

If 

,59 

84-5 

11 

128 

Aerosol 

Channel  1  Data 

58-5,59-1 

12 

129 

II 

tt  2  ** 

59-2,3 

12 

130 

11 

tt  ^  tt 

59-4,5 

12 

131 

tl 

II  ^  It 

60-1,2 

12 

132 

M 

It  ^  It 

60-3,4 

12 

133 

tl 

"  6 

60-5,61-1 

12 

134 

It 

61-2,3 

12 

135 

t; 

"  8 

61-4,5 

12 

136 

tl 

II  g  It 

62-1,2 

12 

137 

tl 

II  JQ 

62-3,4  ' 

12 

138 

II 

"  11 

62-5,63-1 

12 

139 

It 

1.  ^2 

63-2,3 

12 

140 

It 

"  13 

63-4,5 

12 

I 


j 
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A)  The  following  terminology  has  been  used  in  producing  this  document; 

a)  Data  Bank  -  OPAQUE  Random  Access  Disk  Files. 

b)  Data  Base  -  OPAQUE  Data  Managed  by  S2K  DBMS. 

B)  Hyphenated  data  bank  word  numbers  reflect  digit  positions  within  that 
word,  with  position  "1"  always  the  most  significant  digit. 

C)  All  data  base  component  names  have  been  chosen  to  be  Fortran  compatible 
(length  of  name,  default  variable  type)  for  ease  of  use  with  PLI 
software . 


D)  Data  base  component  names  follow  one  of  two  formats  as  follows: 

a)  Non-quality  data 

i)  Letter  "I"  (or  "N"  for  "no.  meas."  data  item)  - 
indicating  Fortran  Integer  type  variable. 

il)  Data  item  code  -  it  should  be  noted  that  when  partial 
data  bank  words  are  used,  an  attempt  has  been  made, 
in  general,  to  indicate  this  by  including  the  digit 
positions  in  the  component  name. 

b)  Quality  data 

i)  Letter  "I"  -  indicating  Fortran  Integer  type  variable, 
li)  Letter  "Q"  -  indicating  quality  data. 

Hi)  Data  base  component  number  -  indicating  which  data 
quality  figure  is  for. 

E)  Data  bank  word  64-2  (data  item  F-2)  is  not  used  and  is  set  to  zero 
in  the  data  bank  -  this  has  no  bearing  on  the  OPAQUE  data  base. 
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DATA  BANK  DATA  BANK  DATA  BASE  DATA  DATA  PICTURE  DECODING 

RECORD  NO.  COMPONENT  NO.  COMPONENT  NAME  ITEM  DESCRIPTION  ALGORITHM 
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FOOTNOTES : 


1)  (Ij^l2l3))  *  (10.  **  Float  (I^  +  6)) 

2)  Iil2^3^4  “  (IJI2I3))  *  (10.  **  Float  (I^  +  8)) 

3)  See  Appendix  1  for  codes  and  their  significance. 

4)  ^2^^2  ~  3)  =  Station  location  code  (see  appendix  2) 

b)  I2  “  Issue  level  of  tape 

i)  1  - 12  ~  Preliminary  data 

li)  ^^l2  *  Reviewed  data 

5)  See  Appendix  3  for  code 

6)  See  Appendix  4  for  code 

7)  Scattering  Code:  =  1  -  Denotes  scattering  being  presented 

=  2  -  Denotes  extinction  being  presented 

*=  3  -  Denotes  both  being  presented 

8)  Filter  Code:  =  4  -  3.4  -  5.5  m 

=  5  -  Open  (2-14  m) 

=6-4  calibration  filter 
=  7  -  To  be  decided 


9) 

Humidity  Code: 

=  1  -  Denotes  direct  humidity  value  to  .1% 

10) 

See  Appendix  4  for 

code 

11) 

See  Appendix  5  for 

code 

12) 

Maximum  Channel: 

=  0  -  No  data 

=  (Max.  channel  -  4)  -  Data  available  (N.B.  - 
data  is  present,  at  least  5  channels  are 
used) . 

13)  ^1^2  Float  (13^12)  (S®®  Appendix  7  for  further  information) 
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14) 

See  Appendix  7  for  code 

15) 

Ijf  8  .125  *  Float  (1^)  ;  I3 

=  9 

No 

measurement 

16) 

13^X2  Integer  value  indicating 

the  nearest 

10°  value  from  true 

North  as  follows: 

-5°  to  5°  -V36 
5°  to  15°  01 


345°  to  355°-^  35 
Special  values: 


00  Calm 

99  No  observa¬ 
tion  or 
failure 


17)  Integer  value  indicating  wind  speed  in  M/sec. 

Ppeclal  value:  ^2^2  '  ^°  failure. 

18)  Ij^l2l2  Integer  value  of  pressure  in  MBAR. 

Ijl2^3  “ 

1^1213  •=  973  .973  MBAR 

(N.B.  Pressure  range  is  approximately  950-1050) 


19) 

^1^2^ 

hV3 

500 

0.1  *  Float 

(I1I2I3) 

'l'2'3 

> 

500 

-0.1  *  Float 

(I3I2I3  -  500) 

h'2'3 

999 

Failure  code 

20)  Iil2’^3  ^^1^2^3^ 

21)  1^1213  0-1  *  float  (I1I2I3) 

Mean  rain  rate  value  in  MM/HR.  The  recommended  sampling 
period  is  4  minutes.  Alternative  periods  must  be  noted 
In  the  comments. 

Special  value:  ~  Failure  code 

22)  See  Appendix  8  for  codes 

23)  I  I_I  Total  rain  in  the  past  hr.  in  MM 

’  Special  value:  1312^3  ^°  failure. 
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APPENDIX  1  -  GROUND  STATE  CODES 

0  =  Dry 

1  =  Wet 

2  =  Flooded 

3  =  Frozen 

A  *  Snow  patches 

5  =  Snow  cover  or  hall  cover 

6  =  Ice  or  frozen  snow  cover 

7  =  Melting  snow 

8  =  Situation  outside  OPAQUE  scheme;  see  comment  list 

9  =  No  observation 


APPENDIX  2  -  STATION  CODES 

1  =  Canada/Denmark 

2  =  France 

3  =  Germany 

A  =  Italy 

5  =  Netherlands 

6  =  UK 

7  =  USA/Gerraany 
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APPENDIX  3  -  UPDATE  STATUS 

NOTE  -  This  Item  does  not  appear  In  the  original  OPAQUE  data  bank,  but 
has  been  included  in  the  data  base  for  convenience.  It  is  meant 
to  serve  as  an  Indicator  of  the  status  of  any  given  logical 
entry  (hour)  in  the  data  base.  The  following  convention  will  be 
used : 

0  =  Raw  data  from  random  access  data  file 

1  = 

2  = 

3  = 

4  = 

5  = 

6  = 

7  = 

8  = 

9  = 
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APPENDIX  _  aerosol  INSTRUMENT  CODE 


0  =  No  data  collected 

1  =  Active  probe  (A) 

2  =  Classical  probe  (C) 

3  =  (A)  and  (C) 

4  =  Royco  (R) 

5  =  (A)  and  (R) 

6  ■=  (C)  and  (R) 

7  =  (A)  and  (C)  and  (R) 

8  =  Other  -  see  comments 


APPENDIX  5  -  AEROSOL  DATA  CODE 

0  *  Data  not  reviewed  for  features 

1  =  Reasonable  approximation 

2  =  Statistically  unreliable 

3  =  Statistically  unreliable  at  larger  sizes  -  low  counts 

4  =  Single  peak  distribution 

5  =  Multiple  peak  exists 
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APPENDIX  6  -  AEROSOL  DATA  DERIVATION 


The  two  digits  of  stored  aerosol  data  are  determined  by  the 
algorithm: 


I  I  -  IFIX(10.*(log( - - -  (N))  +  A)) 

d  log  5 

where  the  following  constraint  applies: 

0  ^  (log( - - -  (N)  +  A)  f  100 

d  log  6 

and  where:  6  =•  Particle  diameter 

N  =  Particle  count/BIN 


In  13  channels  as  Indicated  below: 


0.200 
0.283 
O.AOO 
0.566 
0.800 
1.131 
1.600 
2.263 
3.200 
A. 525 
6.400 
9.051 
12.800 


0.283 
O.AOO 
0.566 
0.800 
1.131 
1.600 
2.263 
3.200 
A. 525 
6.A00 
9.051 
12.800 
18.100 
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APPENDIX  7  -  CONTEL  AND  ADDITIONAX 
MEASUREMENT  CODE 


No  measurement 

Contel 

Sky  camera 

Contel  and  sky  camera 

Additional  measurements  only  as  In  comments 

Contel  and  additional  me.isurements 

Sky  camera  and  additional  measurements 

Contel  and  sky  camera  and  additional  measurements. 
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APPENDIX  8  -  RELIABILIIT  CODES 


First  group,  no  data  available 

0  No  measurement,  as  the  Instrument  is  not  (yet)  working  (min.  period  a 
week,  for  shorter  time  use  1). 

1  No  measurement  in  this  measurement  series  or  may  be  in  the  neighboring 
series  (max.  some  weeks,  for  longer  time  use  0) 
or 

signal  outside  the  range  covered  by  the  instrument  or  signal  destroyed 
by  other  influence  (e.g.  by  looking  direct  in  the  direction  of  the  sun). 


Second  group,  the  relative  error  of  the  data  is  below  50%  (Barnes  15% 
abs.),  i.e.  not  good  for  relationships  but  useful  for  statistics. 

2  Quality  below  normal  standards  is  assumed,  as  the  data  has  not  been 
ciiecked . 

3  Quality  below  normal  standard,  caused  by  measurement  circumstances  (as 
sun  facing  the  instrument  or  fluctuating  light)  or  known  by  detailed 
knowledge  of  the  specific  instrument  (e.g.  noise  at  low  signal  levels), 
or  proven  by  bad  calibration  results. 


Third  group,  the  relative  error  of  the  data  is  below  15%  (Barnes  5% 
abs.),  i.e.  useful  for  all  normal  purposes.  All  data  are  checked  for 
bad  circumstances  as  e.g.  sun  problems. 

4  Quality  normal,  assumed  by  detailed  knowledge,  even  though  it  is  some 
time  since  calibrations  have  been  perfoirmed  or  proven  by  frequent, 
agreeing  calibrations. 
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Fourth  group,  the  relative  error  of  the  data  is  below  5%  (Barnes  2%  abs.)f 
i.e.  they  are  especially  useful  for  demanding  functional  correlation 
relationships . 

5  Quality  above  normal  standard,  based  on  knowledge  of  Instrumental 
performance  and  Information  contained  In  the  basic  data  set  (e.g.  absence 
of  veiling  glare  In  the  path  lumlnence  meters  from  sun  close  to  optical 
axis,  as  derived  from  the  sun  position  and  a  threshold  angle  equivalent 
to  a  5%  effect  with  respect  to  effects  of  artificial  light  sources 

(see  6). 

6  Quality  above  normal  standard,  as  5,  but  also  supported,  for  specific 
interpretations,  by  Information  in  addition  to  that  contained  within 
the  basic  data  set: 

a)  in  case  of  Epply  measurements:  without  any  clouds 
obscuring  the  sun, 

b)  in  case  of  Barnes  measurements:  during  rain  whether  indicated 
or  not  by  rain  rate  meters, 

c)  in  case  of  extinction  coefficient  measurements:  with  values 
at  homogeneous  conditions, 

d)  in  case  of  path  luminance  during  night:  without  any 
occasional  artificial  light  sources  affecting  the  measurement. 
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REDUCTION  AND  ANALYSIS  OF  DATA  FROM  THE  DEFENSE  METEOROLOGICAL  SATELLITE 

PROGRAM 


The  DMSP  program  Is  the  successor  to  the  Air  Force  Defense  Systems 
Application  program  for  performing  meteorological  satellite  imagery  and  to 
the  Data  Acquisition  and  Processing  Program.  DMSP  Is  a  total  program  Involv¬ 
ing  sensors,  data  communications  and  ground  processing  equipment.  Some 
awareness  of  all  of  these  systems  is  necessary  to  understand  and  apply  DMSP 
data.  The  purposes  of  the  DMSP  are  threefold; 

1)  provide  globally  recorded  visual  and  infrared  cloud  cover  and  other 
specialized  environmental  data, 

2)  provide  real  time  direct  readout  of  local  area  environmental  data 
to  mobile  receiving  terminals  throughtout  the  world,  and 

3)  continue  advancement  of  environmental  satellite  technology. 

Basic  interpretation  techniques  are  valid  for  both  conventional  satellite 
data  and  DMSP  data.  But,  DMSP  data  have  several  unique  characteristics  in 
addition  to  a  very  sharp  resolution  capability. 

The  primary  sensors  used  to  obtain  DMSP  data  are  contained  in  an  aero¬ 
space  vehicle  electronics  package  which  is  mounted  on  the  spacecraft  so 
that,  with  nominal  attitude  control,  the  sensors  are  always  oriented  toward 
the  earth.  These  sensors  are  scanning  radiometers  which  respond  to  amounts 
of  radiation  within  specific  spectral  ranges.  Two  types  of  data  are  obtained 
visijal  data,  which  is  a  measure  of  reflected  solar  radiation;  and  Infrared 
data,  which  is  a  measure  of  emitted  earth  and  cloud  radiation. 

The  spectral  interval  used  for  obtaining  visual  data  is  0. 4-1.1  micro¬ 
meters,  l.e.  the  sensor  reacts  as  there  was  a  negative  blue  filter  blocking 
most  of  the  response  in  the  blue  end  of  the  range.  This  minimizes  the  blue 
light  backscatter  of  tht  atmosphere.  The  response  at  the  other  end  of  the 
spectral  range  is  a  result  of  the  type  of  sensor  used.  Its  natural  response 
provides  sensing  into  the  near  Infrared,  but  radiation  sensed  is  still 
primarily  reflected  solar  radiation. 
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A  thermister  bolometer  is  used  to  sense  infrared  radiation  In  the  8-13 
micrometer  range.  This  range  is  used  primarily  because  it  contains  the  peak 
radiations  emitted  by  the  earth  and  its  atmosphere.  The  broad  spectral 
range  was  selected  so  that  the  field  of  view  ot  the  detector  could  be 
narrowed  to  obtain  good  spacial  resolution  while  still  maintaining  a  good 
signal  to  noise  ratio  at  the  detector. 

The  scanning  mirror  in  the  high  resolution,  infrared  mode  (HR/MI) 
rotates  at  1/3  the  speed  of  the  VHR/WHR  mirror,  so  only  1/6  the  number  of 
scan  lines  are  produced.  Tlie  HR/KI  scan  rate  produces  a  2  nautical  mile 
resolution  of  visual  data  and  a  2.4  run  resolution  of  infrared  data. 

Obtaining  proper  spacial  resolution  requires  that  the  tTlR  and  WHR  detectors 
he  smaller  than  their  HR  and  Ml  detectors.  Consequently,  the  VHR  data  are 
available  only  during  daytime  when  there  is  sufficient  reflected  solar 
illumination  of  the  earth  scene  to  provide  adequate  S/N  ratio  at  the  detect-- 
or.  For  infrared  data  there  is  little  diiference  between  the  total  amount 
of  radiation  emitted  from  the  earth  during  day  or  night.  However,  a 
greater  S/N  ratio  problem  exists  for  this  data  than  for  other  detectors. 

The  face  on  the  spacecraft  is  earth  oriented  and  therefore  the  sensors  on 
the  right  are  earth  oriented  and  therefore  the  sensors  on  the  right  are 
earth  oriented.  The  top  of  the  spacecraft  is  always  oriented  toward  the 
sun.  In  this  orbit  there  is  a  possibility,  near  the  poleward  portion  of 
each  revolution,  of  direct  sunlight  impunging  on  the  sensors. 

Some  zero  resolution,  special  purpose,  radiometer  measures  the  amount 
of  sunlight  incident  on  tlie  spacecraft  and,  with  the  earth  oriented  sensor 
system,  measures  the  solar  Illumination  of  the  earth  scene.  This  output  in 
ronjunction  with  a  programmable  gain  memory  unit,  controls  the  signal  level 
output  for  the  visual  data  are  normalized  on  the  spacecraft  for  solar 
illumination;  this  is  an  extremely  important  and  unique  feature  in  data 
prodiictioi’  f:ap.'tbil  i  t  i  es  of  UMSP. 

There  are  four  supplementary  sensors  which  go  on  the  DMSP  spacecraft. 
Normally  each  package  contains  two  supplementary  sensors,  the  exact 
complement  differing  according  to  the  time  period  the  orbit  is  planned  for. 


the  condition  of  the  sensors  on  older  spacecraft,  etc.  One  of  these, 
supplementary  sensor  J  (SSJ)  Is  an  electron  spectrograph  with  one  fixed 
channel  and  one  stepping  channel.  The  channels  detect  energetic  electrons 
over  ranges  of  energies  associated  with  visible  aurora.  The  fixed  channel 
is  6  KeV  and  the  stepping  channel  cycles  through  eight  energy  thresholds: 

54,  98,  219,  600,  1400,  3540,  8200  and  1970  eV.  The  data  sample  is  taken 
approximately  every  second  and  the  field  of  view  is  3  degrees  by  12  degrees. 
SSJ  has  been  Improved;  the  second  generation  is  called  SSJ/2.  It  consists 
of  a  single  stepping  channel  with  six  energy  ranges  with  nominal  energy 
steps  0.3,  0.68,  1.6,  3.5,  7.9,  and  18  KeV,  The  sampling  rate  is  0.0922 
seconds  per  energy  step  and  the  field  of  view  is  a  30  degree  anti-earth 
cone . 


The  J/3  is  such  an  instrument  which  is  an  electron  spectrometer 
designed  for  SAMSO.  The  instrument  consists  of  two  electrostatic 
analyzers  which  measure  the  spectrum  of  precipitating  electrons  that  are 
responsible  for  the  aurora.  The  data  from  J/3  will  be  used  to  determine 
the  location  of  the  boundaries  of  the  auroral  oval  and  will  be  used  to 
correlate  electron  fluxes  with  auroral  observations.  During  the  design 
of  the  instrument,  a  computer  program  was  written  to  trace  the  particle 
orbits  from  the  entrance  slits  to  the  detectors  in  an  effort  to  determine 
the  response  and  detection  efficiencies  of  the  instrument  as  closely  as 
possible.  Tests  Indicate  the  calibration  is  better  than  was  anticipated. 
An  additional  task  that  the  data  will  be  used  for  is  to  correlate  electron 
precipitation  with  substorm  phenomena  at  synchronous  altitudes  along  the 
same  magnetic  field  line  and  for  predicting  total  electron  content  of  the 
ionosphere . 
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Four  types  ot  data  are  produced  by  the  spacecraft's  primary  sensors; 
however,  It  is  not  possible  for  a  tactical  readout  site  to  receive  all  of 
them.  The  DM.Sl'  itansmitter  ope»'ates  at  512  kilobits  per  second  (kbs)  .  The 
normal  operating  mode  is  for  VHR  transmission  on  the  ascending  portion  of 
an  crhit  a:td  WHR  on  the  descending  orbit.  Both  HR  and  Ml  data  are  trans¬ 
mitted  to  all  receiving  sites.  If  the  spacecraft  is  placed  in  an  early 
tr.ordtng  ascending  (c.trly  irv.cning  descending)  orbit  where  there  is  an 
earth  terminator  problem  in  the  visual  scene,  the  '/HR  data  are  transmitted 
cn  the  ascending  portion  of  the  orbit  near  the  equator  where  there  Is 
sufficient  inagiiry.  When  the  spacecraft  progresses  poleward  far  enough 
so  that  the  edge  of  the  VHR  scan  line  enters  the  terminator,  the  trans¬ 
mission  is  switched  to  WHR. 

a.  Visual  Data  The  two  types  of  visual  data,  VHR  and  HR.  complement 
each  other  for  several  meteorological  applicat 1 ons .  The  VHR  is  designed 

to  examine  weather  systems  in  detail  such  that  not  only  is  the  overall 
synoptic  pattern  readily  Identifiable,  but  the  details  of  the  cloud  patterns 
and  their  structure  are  observed.  Since  HR  data  has  a  lower  data  volume 
rate  than  VHR,  it  is  preferred  for  computer  processing. 

b.  Infrared  Data  The  two  types  of  infrared  data  are  Ml  and  WHR. 

The  infrared  sensors  are  co-located  with  visual  sensors.  Coini..on  optics 
insure  that  the  visual  and  infrared  scene  are  viewed  simultaneously  and 

in  the  same  perspective.  This  arrangement  allows  direct  comparison  of  the 
visual  and  infrared  data,  and  permits  the  meteorologist  to  envision  the 
three  dimensional  aspects  of  the  cloud  scene  generated  by  atmospheric 
motions. 

The  MI  and  WHR  sensors  were  developed  to  serve  several  different 
purposes  for  the  meteorological  analyst.  MI  data,  with  its  2. A  nra  spatial 
resolution  at  subpoint,  is  used  as  the  primary  infrared  data  for  comparison 
with  all  visual  dat.a.  It  is  also  the  more  accurate  infrared  data  for  use 
in  quantlta. ivf  thermal  determinations. 
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D^ta  The  DMSP  sensors  were  designed  to  provide  high  quality  data. 
The  system  approach  was  taken  to  ensure  that  the  resolution  inherent  In 
the  basic  sensor  electronic  signal  would  be  faithfully  preserved  through 
amplification,  recording,  transmission,  and  eventual  display  on  the  ground. 
The  sensors  were  engineered  to  provide  subpoint  spatial  resolution  of  1/3 
NM  for  VHR  and  WHR  data,  2  NM  for  HR  data,  and  2.4  MM  for  MI  data. 

All  these  data  are  rectified  by  the  Signal  Processor  No.  1  but  the 
rectification  process  does  not  improve  spatial  resolution. 

As  a  general  rule  of  thumb,  VHR  imagery  displayed  in  a  1:15  million 
scale  has  a  spatial  resolution  of  one-third  to  one-half  nm  within  the  first 
inch  of  data  either  side  of  subpoint.  From  one  inch  and  two  and  one-half 
inches  either  side  of  the  subpoint,  the  spatial  resolution  degrades  from 
one-half  to  one  nm.  From  the  two  and  one-half  inch  point  to  data  edge, 
the  spatial  resolution  further  degrades  to  two  nm.  These  degradations 
are  due  to  foreshortening. 

Spatial  resolution  is  more  difficult  to  quantitatively  determine  for 
infrared  sensors  than  it  is  for  visual  sensors.  In  the  infrared,  the 
spatial  resolving  power  of  the  sensor  can  be  masked  by  the  thermal 
accuracy  of  the  sensing  system. 

•  The  error  between  actual  temperature  and  effective  radiating  tempera¬ 
ture  was  observed  to  be  related  to  the  moisture  distribution.  The 
regression  equation  describing  it  is: 

Error  (T^  -  ERT)  =  35.97  +  30.98  Log  (w/cos3) 

where:  w/Cos&  represents  a  line  of  sight  measurement  through  the  moisture 

] ayers . 
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A  study  on  absolute  thermal  accuracy  was  performed.  Radiosonde  data 
was  used  from  within  six  hours  of  a  satellite  pass.  A  correlation  was 
derived  between  radiosonde  "observed"  cloud  top  temperatures  and  the 
effective  radiating  temperature  Indicated  by  MI  data. 

It  Is  apparent  that  atmospheric  attenuation  factors  are  different 
for  each  site  depending  on  its  location  and  season.  Many  sites  have  a 
continuing  evaluation  program  so  that  accurate  temperatures  can  be  obtained 
for  mission  support.  To  assist  those  who  do  not  have  access  to  site 
peculiar  calibration  data,  theoretical  atmospheric  attenuation  curves  for 
a  standard  atmosphere  have  been  developed  for  the  MI  sensor. 

ERT,  ,  -  ERT  +  0.0125  (ERT  -  220)  +  0.00156  (ERT  -  220)^ 

(corr) 

Orbit  Characteristics  The  ability  to  use  DMSP  data  Is  greatly 
enhanced  by  unique  data  location  procedures.  Data  location,  or  gridding 
of  data,  at  the  very  best,  is  a  complex  science.  An  understanding  is 
necessary  of  the  sensing  mechanism,  spacecraft  orbit,  and  display  equip¬ 
ment.  The  scanners  are  essentially  looking  at  a  curved  surface  from  a 
fixed  point  in  space.  The  resulting  effect  is  termed  foreshortening.  Not 
only  are  these  data  foreshortened  because  of  "looking"  over  the  curvature 
of  the  earth,  but  the  spacecraft  travels  in  an  orbit  such  that  the  longi¬ 
tude  and  latitude  lines  are  curved  and  their  orientation  changes  depending 
on  the  subpoint  latitude. 

The  nominal  orbit  for  the  DMSP  spacecraft  is  a  450  nm  circular  orbit 
with  an  inclination  angle  of  98.7  (exactly  98.747).  This  means  the  orbital 
plane  Is  inclined  98.7  to  the  equatorial  plane  where  the  spacecraft  crosses 
the  equator  northbound  (ascending  node).  The  inclination  angle  was 
selected  to  Insure  that  the  nominal  450  nm  circular  orbit  is  sun-synchronous 
(maintains  a  relatively  constant  relationship  to  the  sun  so  that  the 
ascending  node  remains  at  a  constant  solar  time).  The  combination  of  a 
450  nm  circular  orbit  and  98.7  inclination  angle  means  that  the  satellite's 
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orbital  plane  rotates  slowly  around  the  earth  at  the  same  rate  and 
direction  that  the  earth  rotates  around  the  sun.  At  the  most  poleward 
positions  of  the  spacecraft,  its  motion  is  tangent  to  the  parallels  of 
latitude,  adding  to  the  complexity  of  a  grid's  latitude/ longitude  configura¬ 
tion. 


The  nodal  period  of  this  sun-synchronous  orbit  is  101.56  minutes. 
During  the  time  of  one  revolution,  the  earth  revolves  under  the  orbit  to 
the  east  and  the  orbital  plane  precesses  slightly  to  the  east.  The  result 
is  that  each  nodal  crossing  is  approximately  25.4  west  of  the  previous 
crossing.  Approximately  every  five  days  there  is  a  15  revolution  day; 
therefore,  a  space-craft  peculiar  cycle  is  established  which  causes  both 
dally  data  coverage  and  local  time  of  data  variations. 

Once  DMSP  data  are  acquired,  only  the  longitude  and  nodal  crossing  and 
the  elapsed  time  since  the  nodal  crossing  are  needed  in  order  to  accurately 
locate  these  data. 

At  AFGWC  a  continuous  sinusoidal  correction  is  applied  to  the  imagery. 
A  tactical  site  usually  sets  an  altitude  setting  which  is  assumed  to  be 
valid  for  the  few  minutes  the  spacecraft  is  within  readout  range  of  the 
site. 


Data  is  received  from  various  satellites  as  they  orbit  the  earth. 

This  data  is  archived  onto  tape  by  the  World  Data  Center  A  for  Solar- 
Terrestrial  Physics  in  Boulder,  Colorado.  For  the  past  eighteen  months 
and  on  a  continuing  basis  they  have  been  and  are  sending  tapes  containing 
the  data  collected  from  two  satellites  F2  and  F4 .  Various  other  satellites 
will  be  sent  up  in  the  future  and  their  data  in  turn  will  be  processed  in 
the  same  manner  as  F2  and  F4,  (though  at  a  higher  data  rate). 

One  magnetic  tape  from  Boulder  contains  approximately  ten  days  worth 
of  data.  If  there  are  two  satellites  (as  at  present)  on  a  tape  each  tape 
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contains  five  calendar  days  worth  of  data.  After  considerable  effort  on  our 
part  to  make  the  existing  processing  programs  as  efficient  and  compact 
as  possible,  it  still  takes  approximately  2000  CP  seconds  computer  time  to 
set  up  the  present  data  base  and  produce  the  microfiche  flux  plotting  for 
five  calendar  days  from  one  satellite.  This  data  base  consists  of  one 
tape  per  15  calendar  days  of  packed  channel  counts  for  one  satellite  and 
one  tape  per  month  of  flux  data.  Already  under  this  present  system 
approximately  two  years  worth  of  data  of  one  satellite  and  six  months  of 
two  satellites  have  been  processed  and  we  now  have  a  data  base  of  approxi¬ 
mately  eighty  tapes.  Some  of  the  tapes  are  labeled  and  some  not.  We  have 
found  several  tapes  mislabeled  and  others  with  data  missing.  Also  this 
processing  will  expand  to  accommodate  three,  four  or  more  satellites  sent 
up  in  the  future. 

The  preceedlng  describes  the  software  required  in  setting  up  a  CDC 
data  base.  This  consists  of  only  a  part  of  the  total  project.  On  a 
continuous  day  to  day  basis,  we  receive  several  requests  to  run  several  of 
the  existing  application  programs  which  keep  a  programmer  quite  busy.  In 
addition  it  involves  the  running  of  the  statistical  survey  program  whose 
input  is  all  of  the  data  base  flux  tapes  (approximately  40  tapes). 

Through  meticulous  effort  put  into  packing  the  data  and  by  use  of 
high  density  tapes,  we  feel  we  can  put  about  six  to  eight  months  worth  of 
data  onto  one  tape.  The  old  system  would  require  18  to  24  tapes.  The 
production  runs  will  be  reduced  to  a  fraction  of  what  they  are  now  (less 
than  500  seconds).  We  intend  to  set  up  a  partly  interactive  system  whereby 
an  operator  at  an  Intercom-terminal  will  type  in  the  date  of  data  desired, 
the  time  Intervals,  and  application  program  desired  and  a  batch  job  would 
be  generated  to  run  the  request.  After  logging  in  the  operator  would  be 
prompted  by  clear  and  straightforward  messages.  Any  of  the  countless  day- 
to-day  requests  could  be  generated  by  the  scientists  or  Initiators  them¬ 
selves  even  if  they  have  little  or  no  programming  experience.  We  would 
then  have  more  time  to  develop  new  software  which  could  be  very  easily 
added  into  the  existing  software. 
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REDUCTION  AND  ANALYSIS  OF  HIGH  ENERGY  ELECTRON  DATA 
DMSP  PROGRAM 


1.  Normal  production  of  microfiche  and 
archived  output  tapes. 

2.  Production  of  particle  spectrum  and  paper 
plots  for  publication. 

3.  Dump  of  the  edge  of  the  auroral  oval. 

U.  Plar  maps  of  auroral  oval. 

5.  Replottlng  microfiche  with  abnormal 
problems . 

6.  Testing  method  for  evaluating  equatorial 
edge  of  oval. 

7.  Checking  particle  data  against  data  from 
existing  DMSP  tapes. 
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SUMMARY 


1.  Tape  Output 

Files  are  accumulated  on  tapes  up  to  six  production  runs 
on  one  output  tape.  Thus,  one  production  run  can  affect 
others.  Sometimes  Jobs  abort  through  mis-keypunching , 
parity,  tape  breaking,  system  problems,  etc.  and  thus 
garbage  data  and  half  files  are  put  on  tape.  The  further 
runs  have  to  read  data  on  tape  before  adding  data  and  runs 
are  not  independent.  The  ideal  is  to  get  consistently 
one  run  a  night  for  each  calendar  months  processing.  Large 
chunks  of  time  of  say  5,  6  or  8  hours  are  not  too  beneficial 
as  it  now  stands. 

2.  Disk  Output  -  Achieved  to  Tape 

This  method  is  more  appropriate  if  we  are  getting  discrete 
chunks  of  say  3,  6,  or  8  hours  of  computer  time.  Runs  are 
completely  independent.  Output  files  are  catalogued  on 
disk.  However  this  Involves  the  extra  Step  of  achieving 
disk  files  to  tape  then  purging  the  disk  file.  About  four 
production  runs  could  be  put  on  a  disk  pack.  To  make  this 
feasible  we  would  need  another  disk  pack. 
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Typical  DMSP  Activities 


Generate  high  latitude  maps  of  the  electron  precipitation 
patterns  based  upon  Kp. 

Using  the  SSF/3  data  determine  the  relationship  between  Kp, 
the  equatorward  boundary  of  the  auroral  oval,  and  magnetic 
local  time.  This  requires  the  determination  of  the 
equatorward  boundary  for  several  thousand  orbits  of  SSJ/3 
data. 

The  project  is  attempting  to  find  a  way  to  rapidly 
determine  the  equatorward  edge  of  the  oval  and  to  specify 
Its  position.  The  work  already  done  on  the  large  statisti¬ 
cal  program  has  shown  that  a  fully  automated  determination 
of  the  boundary  is  not  possible.  Hence,  some  interactive 
capability  perhaps  using  the  Tektronix  graphics  is 
required. 

Determine  an  algorithm  to  be  used  in  real  time  at  AWS  for 
giving  the  approximate  position  of  the  equatorward  boundary. 
This  project  requires  the  rtudy  of  a  number  of  auroral 
passes  of  SSJ/3  and  come  up  with  several  possible  algorithm. 

Over  the  long  run  AWS  would  like  to  see  the  global 
statistical  study  redone  based  on  different  energy  bins. 

Now,  in  the  global  study  the  integral  flux,  energy  flux 
and  average  energy  are  the  three  quantities  studied.  For 
the  follow  on  study,  they  would  like  to  see  the  same 
quantities  studied  but  to  have  them  calculated  separately 
for  the  energy  ranges  causing  the  F  &  E  regions  in  the 
ionosphere  respectively.  In  addition,  they  would  like  to 
see  the  global  maps  combined  with  a  simple  code  to  actually 
predict  the  E  region  profile. 
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Finally,  AWS  would  like  AFCL  to  produce  an  auroral  activity 
index  using  the  equatorward  boundary  of  the  oval.  This 
would  require  that  for  every  pass  of  the  Si»J/3  data  the 
boundary  be  determined. 

Many  of  these  tasks  have  been  completed.  The  statistical  program  is 
finished.  The  development  of  the  plotting  programs  to  display  the  data 
must  still  be  completed. 

The  algorithm  to  find  the  equatorward  boundary  has  also  been  success¬ 
fully  tested.  This  work  is  being  done  to  support  the  Air  Weather  Service. 
Since  Air  Weather  Service  plans  to  use  the  algorithm,  the  programming  had 
to  be  fully  documented.  This  will  be  Incorporated  in  an  in-house  report. 

The  printout  program  is  also  working.  The  test  we  use  to  determine 
where  to  printout  the  data  is  not  as  accurate  as  the  algorithm  for  determin¬ 
ing  the  equatorward  boundary.  This  program  was  modified  to  use  the 
algorithm  to  determine  the  equatorward  edge  and  to  printout  the  data  for  1 
minute  before  and  after  this  points 

Additional  work  falls  into  two  categories.  The  first  category 
includes  projects  to  modify  and  extend  the  global  statistics  work.  The 
second  category  includes  projects  to  analyze  the  boundary  determinations  we 
have  made  using  the  SSJ/3  data. 

The  first  project  in  extending  the  global  statistics  work  is  to  change 
the  binning  size.  At  present  the  bins  are  one  half  hour  wide  in  magnetic 
local  time  by  2  degrees  in  magnetic  latitude.  This  means  that  as  a  function 
of  increasing  latitude,  the  area  of  the  bins  decrease  and  the  total  number 
of  points  averaged  In  each  bin  is  not  the  same.  To  partially  eliminate 
this  problem,  the  bins  should  vary  in  latitudinal  size  as  follows; 
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50-60  degrees 
60-70  degrees 
70-80  degrees 
80-90  degrees 


1/2  degree  magnetic  latitude  x  1/2  hr  MLT 

1  degree  magnetic  latitude  x  1/2  hr  MLT 

2  degrees  magnetic  latitude  x  1/2  hr  MLT 
2.5  degrees  magnetic  latitude  x  1/2  hr  MLT 


The  averaging  Interval  should  also  be  decreased  from  10  seconds  to  4 
seconds.  This  will  require  generating  a  tape  of  4  second  averages.  The 
program  to  Ignore  bad  data  will  have  to  be  modified  to  run  off  such  a  tape. 
We  will  wish  to  redo  the  statistical  study  using  two  years  of  data  with  the 
new  binning  and  the  new  averaging  intervals. 


The  second  new  project 'involving  global  statistics  involves  binning 
the  data  according  to  the  orientation  of  the  Interplanetary  Magnetic  Field. 
We  have  6  months  of  data  on  the  IMF.  The  data  exist  as  1  hour  averages 
available  on  tape  and  listing.  Because  of  the  crude  time  averages  of  the 
IMF  we  have,  only  the  original  average  sizes  and  bin  sizes  need  to  be  used. 
The  programming  will  Involve  the  capability  to  determine  for  each  ten 
second  average  the  orientation  of  the  IMF  for  a  period  1  hour  before  the 
time  of  the  ten  second  average.  As  a  first  division  the  data  will  be 
binned  as  follows: 
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Finer  division  according  to  the  actual  magnitude  of  B  and  B  will  be  an 

z  y 

additional  step  after  this  initial  work  is  completed. 
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As  a  parallel  job  wo  have  been  determining  the  equatorward  boundary 
of  the  auroral  zone  from  the  survey  microfiche.  Approximately  5000  to 
6000  boundaries  have  been  calculated  and  typed  onto  computer  cards.  Kach 
card  includes  the  boundaries  on  each  side  of  the  oval  for  each  pass  over 
the  north  or  south  pole  (magnetic  latitude  and  magnetic  local  time).  They 
include  also  the  day  and  year  of  each  crossing,  the  approximate  time  in 
the  clay  in  Universal  Time  seconds,  and  the  Kp  value  for  the  three  hour 
interval  in  which  measurement  was  made. 

We  anticipated  first  analyzing  the  manner  in  which  the  equatorward 
boundary  varies  in  magnetic  latitude  and  local  time  as  a  function  of  the 
geomagnetic  activity  index  Kp.  Specifically,  we  dccslgned  programs  with 
the  capability  to  read  each  card,  determine  in  what  hour  of  magnetic  local 
time  the  boundaries  occur  and  for  each  whole  hour  in  magnetic  local  time 
to  bin  the  magnetic  latitude  of  the  boundary  according  to  each  level  in  Kp. 
We  have  calculated  approximately  twnety  boundaries  for  each  level  of  Kp  (0, 

0+,  1- . etc.)  in  each  hour  of  local  time  to  which  the  satellite 

orvit  gave  access.  For  each  level  of  Kp  in  a  magnetic  local  time  zone  the 
program  calculates  the  average  and  standard  deviation  of  the  magnetic 
latitude  of  the  boundaries.  The  program  also  calculates  using  all  of  the 
boundary  determinations  in  a  local  time  zone,  the  least  squares  straight 
line  fit  to  the  data  and  the  linear  correlation  coefficient. 

Secondly,  we  anticipated  analyzing  the  manner  in  which  the  morning 
equatorward  boundary  varies  as  a  function  of  the  evening  equatorward 
boundary.  Each  computer  card  has  an  evening  and  morning  equatorward 
value  for  a  pass  over  the  north  or  south  pole.  What  we  require  are  programs 
that  will  read  the  cards  and,  for  specific  magnetic  local  time  bins  on  both 
the  evening  and  morning  sides  of  the  oval,  will  plot  the  magnetic  latitude 
of  the  evening  side  boundary  versus  the  magnetic  latitude  of  the  morning 
side  boundary.  For  example,  we  may  need  to  look  at  those  orbits  that  pass 
through  the  magnetic  local  time  zone  1800  to  1900  on  the  evening  side  and 
0500  to  0600  on  the  morning  side.  For  each  orbit,  so  Identified,  we  wish 
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to  plot  the  latitude  of  evening  boundary  vs  the  latitude  of  the  morning 
boundary.  For  each  pair  of  local  times  selected  all  the  points  should  be 
plotted  and  the  data  least  square  fitted  and  the  linear  correlation  co¬ 
efficient  calculated. 

The  last  project  is  a  study  of  the  manner  in  which  the  equatorward 
boundary  varies  as  a  function  of  the  orientation  of  the  IMF.  As  in  the 
other  statistical  study  using  the  IMF,  this  project  requires  the  ability 
to  read  the  cards  of  the  IMF.  For  each  one  hour  magnetic  local  time  zone, 
the  program  must  determine  the  manner  in  which  the  equatorward  boundary 
varies  as  a  function  of  the  and  components  of  the  IMF  according  to 
the  following  separations. 

a. 

7  0. 

b. 

B  ✓  0. 
z  ^ 

c . 

all  B  . 

y 

For  each  bin  in  B  or  B  ,  the  mean  and  standard  deviation  should  be 
calculated  and  using  all  the  data  points  the  least  square  fit  and  linear 
correlation  coefficient  should  be  calculated.  As  in  the  other  study  with 
IMF,  the  1  hour  average  should  be  used  that  proceeds  by  1  hour  the  time 
at  which  a  given  boundary  determination  was  made.  The  study  also  might  be 
later  expanded  to  include  separations  according  to  the  solar  wind  plasma 
velocity . 


B^  in  1  gamma  bins  from  +15  gamma  to  -15  gamma  for 
B^  in  1  gamma  bins  from  +15  gamma  to  -15  gamma  for 
B^  in  1  gamma  bins  from  +10  gamma  to  -10  gamma  for 
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DMSP  J/3  Processing:  Additional  tasks 

Priority  1  Priority  2 

Electron  Precipitation  high 

latitude  maps:  algorithm: 


Contour 

Plot 


Contour  plots 
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J/3  Processing:  Additional  Tasks  (continued) 


Priority  3 

Auroral  Oval  Determination 


DMSP  Stacked  Data  Tape  (DSDT)  Format 


1.  General  Format 

The  DSDT  tape  contains  a  distillation  of  the  original  data  tapes 
written  by  the  Digital  Data  Formatting  System  (DDFS) .  Bad  data  has  been 
eliminated,  and  the  data  has  been  separated  out  and  calibrated. 

The  DSDT  tapes  are  written  at  556  BPI,  in  odd  parity,  on  a  7  track 
drive. 

The  DSDT  tapes  contain  a  series  of  records  making  up  on  giant  file 
containing  all  of  the  data  from  a  single  pass.  (This  file  may  span  2 
tapes).  Each  record  contains  one  data  line.  Three  types  of  records 
(data  lines)  appear:  V  Data,  I  Data,  and  T  Data.  T  Data  being  calibrated 
I  Data.  Generally,  but  not  necessarily,  due  to  noise  problems  in  the 
original  data  which  may  cause  loss  of  one  or  more  data  lines,  6  V  records 
will  appear  for  each  I  and  T  record  pair,  {k  T  record  always  follows  each 
I  record). 

Each  record  consists  of  an  Identifier  portion  and  a  data  portion.  The 
identifier  portion  contains  twnety  identifier  parameters.  The  first 
parameter  is  a  record  (data)  type  discriminator,  4  BCD  characters,  WW  for 
V  Data,  nil  for  1  Data,  and  TTTT  for  T  Data.  (Each  BCD  character  -  6  BITS). 

The  remaining  nineteen  identifier  parameters  are  24  BIT  binary 
integers  (4  tape  characters  each). 

The  data  values  are  either  6  BIT  digital  count  values  (V  Data  and  I 
Data)  or  24  BIT  Integers  (T  Data)  as  described  below. 


II.  V  Data 


T  =  63 
I  =«  31 
V  =  65 
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Identifier  Parameters 


1.  WW 

2.  V  line  number  -  counts  V  lines  written  on  DSDT  tape. 

3.  Pass  number 

A.  DMSP  satellite  number 

5.  Data  day 

6 .  Data  month 

7.  Data  year 

8.  First  line  time,  tenths  of  seconds  after  ascending  node 

9.  Current  line  time,  tenths  of  seconds  after  ascending  node  (counting 
up  from  first  line  time  described  above). 

10.  DZ  -  zero  resolution  sensor  output  (0-63). 

11.  DR  -  dark  reference  value  (0-63). 

12.  PSGC  -  preset  gain  code  (0-63)  (not  yet  used). 

13.  Bad  sample  count  -  number  of  samples  set  to  zero  due  to  problems  with 
the  data. 

14.  Subpoint  sample  number  (not  yet  used). 

15.  Satellite  heighi.,  N.  miles  (not  yet  used). 

16-20, Not  yet  used. 

Following  the  20th  identifier  parameter  are  3400  V  Data  values,  each 
a  6  BIT  binary  Integer  (tape  character).  A  zero  corresponds  to  black,  63 
to  white.  The  data  includes  the  complete  earth  scan. 
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III.  I  Data 


Identifier  Parameters 

1.  nil 

2.  I  line  number  -  counts  I  lines  written  on  DSDT  tapes. 

3-9.  Same  as  V  Format 

10.  Thermistor  output  count  (0-63). 

11.  Calibration  source  IR  count  (0-63). 

12.  Free  space  reference  count  (0-63). 

13.  31  OK  reference  count  (0-63). 

lA.  21  OK  reference  count  (0-63). 

15.  First  horizon  sample  number. 

16.  Second  horizon  sample  number. 

17.  DM  data  word  -  not  yet  used. 

18.  Satellite  height,  N.  miles  -  not  yet  used. 

19-20.  Not  used. 

Following  the  20th  identifier  parameter  are  700  T  Data  values.  Each 
is  a  6  BIT  binary  integer  (tape  character).  0  represents  the  warmest  value 
and  63  the  coldest.  The  data  includes  the  complete  earth  scan. 

IV.  T  Data 

Identifier  Parameters 

1 .  TTTT 

2-20. Same  as  I  Data 

Following  the  20th  Identifier  parameter  are  700  calibrated  temperature 
values  corresponding  exactly  to  the  700  I  Data  values  in  the  preceedlng 
I  Data  record.  Each  temperature  value  is  a  2A  BIT  binary  integer  (A  tape 
characters).  The  values  are  in  tenths  of  a  degree  celslus.  Negative 
values  follow  the  convention  by  which  77777776g  is  1,  77777775g  is  2, 
etc. ,  and  00000000.  or  77777777.  is  0. 

O  O 
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DATA  FORMATTER  OUTPUT  FORMATS  FOR  DMSP  BLOCK  5C  PRIMARY  DATA 

I.  FORMATS 

1.  Data  Is  arranged  In  36  binary  bit  words  on  seven  track  tapes  with 
odd  parity  and  have  a  density  of  800  bits  per  inch. 

2.  The  Data  h'orraatter  output,  or  (DF)  data  tapes  are  generated  in 
an  IBM-compatible  format.  Each  DF  tape  contains  a  header  record,  N  number 
of  data  records,  and  an  End-of-Flle,  followed  by  another  header  record. 

When  data  runs  require  two  DF  tapes,  the  first  tape  ends  with  two  End-of- 
Flles,  and  the  second  data  tape  eands  with  an  End-of-File  followed  by  the 
header  record.  The  header  record  contains  27  words  and  provides  unique 
identification  for  each  data  run. 

3.  Block  30  primary  data  can  be  in  either  of  two  different  formats: 
Analog  or  Digital.  Column  7  of  the  header  record  contains  the  control  that 
designates  either  Analog  or  Digital  data:  An  "A"  in  this  column  designates 
Analog:  a  Blank  designates  Digital. 

a.  In  the  "Analog  Mode",  an  analog  sub-carrier  is  aemodula- 
ted  and  digitized  by  the  DF; 

b.  In  the  "Digital  Mode",  the  analog  subcarrier  is 
processed  by  a  Block  V  Formatter  (BVF)  at  the  Command 
Readout  Station  (CRS)  and  sent  to  AFGWC  via  an 
American  Satellite  Communications  link  as  Digital  Data, 
where  the  DF  then  formats  the  data  for  processing. 

The  analog  video  base  bands  may  be  digitized  by  either 
the  DF  in  Analog  Mode  or  a  BVF  in  Digital  Mode. 

II.  ANALOG  MODE 


1,  In  the  Analog  Mode,  Il/I  data  is  formatted  into  38^0  (decimal) 
word  records  for  processing.  The  output  formats  for  this  Analog  Mode  data 
are  shown  in  Figures  5  and  6.  It  should  be  noted  that  the  3840  word  record 
length  provides  for  an  Integral  number  of  H  and  I  scan  lines  in  each  record 
(no  partial  scan  lines).  The  first  word  found  in  each  record  will  be  a 
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"frame/line  counter"  and  the  last  word  will  always  be  an  "assembled  time- 
word"  . 


2.  Each  scan  line  Is  provided  with  a  unique  frame/ line  counter  to 
designate  the  beginning  of  the  scan  line.  As  part  of  this  frame/line 
counter,  the  most  significant  bit  provides  an  indicator  flag  for  the  type 
of  data  to  follow.  A  '1'  in  this  position  indicates  visible  data  (H  or  V) : 
a  'O'  indicates  infrared  data  (I  or  W) . 


FRAME/LINE 

CALIBRATION 

ii/i 

CALIBRATION 

TIME 

COL'NTEP 

DATA 

DATA 

DATA 

BLANKS 

WORDS 

"ANALOG"  SCAN  LINE  FORMAT 

3.  The  calibration  slots  are  digitized  and  depicted  by  several 
calibration  words.  These  calibration  words  are  repeats  of  the  data  in  the 
calibration  slots;  thus,  if  the  time  code  slot  contained  a  '1'  bit,  the 
time  code  words  would  contain  77  (octal).  In  the  "Digital  Mode"  these 
calibration  slots  have  been  assembled  into  a  single  calibration  word  for 
the  several  calibration  parameters. 

4.  The  H  and  I  data  may  be  in  line  pairs;  that  is,  a  visible  scan 
line  followed  by  an  infrared  scan  line.  Or  the  data  may  be  H  only  or  I 
only,  but  the  same  record  length  is  still  used.  The  most  significant  bit 
in  the  frame/line  counter  is  the  key  to  identifying  the  data. 

5.  At  the  end  of  each  scan  line  there  is  an  "assembled  spacecraft 
time-word".  This  time-word  is  updated  onboard  the  spacecraft  once  every 
two  minutes  real-time.  Therefore,  this  time-word  (in  the  H  and  I  data) 
will  change  approximately  once  every  214  scan  lines.  However,  if  an  error 
occurs  during  the  assembly  of  a  time-word,  the  previous  value  will  be  held 
and  the  erroneous  time-word  discarded  until  the  next  two-minute  update 
occurs.  Thus,  it  is  possible  to  have  more  than  214  scan  lines  between 
time-words  (  a  multiple  of  214),  however,  in  most  cases,  there  will  be 
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about  214  scan  lines.  In  the  V/W  mode,  the  time-word  changes  approximately 
every  1284  scan  lines. 

6.  The  V/W  data  is  similar  to  the  H/1,  except  the  record  length  is 
3222  words  (decimal),  and  the  data  is  collected  so  tViat  blocks  of  data  are 
processed  on  the  tape.  Due  to  the  high  density  of  the  V/W  data  and 
transfer  rate  constraints  of  the  tape  control  unit,  V/W  data  is  processed  in 
one  of  two  different  modes:  Only  every  other  scan  line  is  processed;  or 
only  one  scan  line  In  four  is  processed  (as  designated  by  the  DF  set-up 
card,  Column  10. 

It  is  possible  for  a  switch  from  visible  to  IR  data,  or  IR  to  visible 
data,  to  occur  within  any  given  data  run  (watch  for  the  significant  bit  in 
the  Frame/Line  Counter). 

III.  DIGITAl.  MODE 


1.  The  digital  data  is  processed  in  a  similar  manner,  except  the 
calibration  information  has  been  placed  in  one  word  at  the  end  of  the  data 
stream.  Thus,  the  total  length  of  each  scan  line  has  been  reduced  by 
removing  the  calibration  overhead.  The  exact  record  lengths  for  processing 
the  digital  H/I  and  V/W  data  have  not  been  selected;  however,  if  a  digital 
tape  header  record  is  examined.  Columns  2  through  6  contain  the  number  of 
words  in  each  record  will  be  adjusted  to  provide  an  integral  number  of 
scan  lines  per  record. 


FRAME/LINE 

CALIBRATION 

TIME 

COUNTER 

DATA 

WORD 

WORD 

"DIGITAL"  SCAN  LINE  FORMAT 
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IV.  THE  DATA 


1.  An  Individual  data  sample  consists  of  raw  radiances  as  detected 
by  the  onboar  sensor  In  six-bit  binary  representation  of  the  radiance 
value,  where  zero  is  dark  and  77  (octal)  is  maximum  brightness  in  visible 
data.  The  IR  data  sample  is  similar;  however,  may  seem  reversed  in  a 
depiction.  The  zero  is  210  degrees  Kelvin  and  77  (octal  is  310  degrees 
Kelvin. 

V.  UNIVAC  TAPE  FORMATS 


1.  All  records  are  separated  by  an  inner  record  gap  referred  to  as 
an  End-Of-Record  (EOR)  mark.  The  EOR  is  a  s/4  inch  blank  gap  on  the  tape. 
When  an  I/O  request  is  made  to  a  tape  unit,  the  entire  data  record  (EOR  to 
EOR)  is  read/written  by  the  tape  unit.  If  the  user  requests  only  1  word  of 
a  2000  word  record,  1  word  will  be  transferred  to  his  input  buffer  but  the 
tape  will  advance  to  the  beginning  of  the  next  record,  i.e.  the  remaining 
1999  words  of  the  record  will  be  dropped. 

2.  All  files  on  a  tape  are  separated  by  End-Of-File  (EOF)  marks. 

The  EOF  is  a  hardware  mark  written  on  the  tape  by  the  tape  drive.  All  tape 
displays  use  at  a  minimum  2  files.  The  first  file  contains  descriptive 
information  and  data  records.  The  second  file  usually  contains  an  End-Of- 
Display  sentinel.  On  occasion,  it  may  contain  a  continuation  sentinel 
which  indicates  that  the  display  is  continued  on  another  tape.  The  data 
file  contains  integral  numbers  of  scan  lines  which  are  used  to  build  a 
picture  a  line  at  a  time.  The  number  of  scan  lines  per  record  is  constant 
for  any  given  display,  but  may  vary  from  display  to  display  because  of  the 
different  map  scales  employed. 

VI.  DATA  FORMATTER  TAPES 


Tapes  produced  by  the  Harris  Corporation  Data  Formatter  (DF)  contain 
the  following  information: 


1.  Record  1  (DF  Header)  The  first  record  contains  28  words  in 
column  binary  format.  A  breakdown  of  column  binary  is  given  in 
attachment  1.  The  8  column  DF  card  input  from  which  the  28  word  record  is 
derived  is  given  In  attachment  2.  This  card  Is  read  by  the  DF  In  column 
binary  and  written  to  record  number  1  of  the  tape.  The  processing  para¬ 
meters  specified  are  used  to  direct  both  the  data  formatter  and  subsequent 
computer  processing  on  a  larger  mainframe.  Columns  2-6  of  the  card  specify 
the  size  of  data  records  desired.  This  entry  is  a  5  digit  octal  number. 

2.  Record  2  through  N  (DATA  NNNNN  WORDS)  All  subsequent  records 
on  the  DF  tape  contain  NNNNN  data  words  where  NNNNN  is  the  record  size  as 
specified  on  the  DF  card.  Each  record  after  the  header  record  up  to  and 
including  the  End-Of-File  mark  (EOF)  contains  a  fixed  number  of  words. 

The  size  of  each  record  is  specified  by  columns  2-6  of  the  DF  card.  In 
addition,  each  record  of  NNNNN  words  contains  an  integral  number  of  scan 
lines.  Video  and  IR  data  lines,  when  both  are  present,  are  interleaved. 

Note  that  video  lines  are  distinguished  from  IR  lines  by  a  frame  line 
counter  which  has  bit  35  (leftmost  bit)  set. 
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